CAPABILITY MEASUREMENTS FOR MULTIVARIATE PROCESSES: DEFINITIONS AND AN EXAMPLE FOR A GEAR CARRIER Steven J. Littig C. Teresa Lam and Stephen M. Pollock Department of Industrial & Operations Engineering University of Michigan Ann Arbor, MI 48109-2117 Technical Report 92-42 June 1992 Submitted to HE Transactions

Capability Measurements for Multivariate Processes Definitions and an Example for a Gear Carrier Steven J. Littig, C. Teresa Lam and Stephen M. Pollock Department of Industrial and Operations Engineering The University of Michigan, Ann Arbor, MI 48109 Abstract A new approach to process capability measurement is developed and applied to a complex multivariate tolerance system. The approach, including a comprehensive and complementary system of indices, is based on process proportion nonconforming and deviation from target. It is both consistent and extendible for any univariate or multivariate tolerance systems. A detailed example for a transmission gear carrier is presented to highlight the solution method for single holes, coaxial holes, and multiple coaxial holes. KEY WORDS: Coaxial hole tolerance; Hole location measurements; Multivariate capability; Multivariate tolerances; Process capability. 1. Introduction Two critical elements needed for continuous improvement in manufacturing industries are the ability to quantitatively measure production quality, and the use of those measures for competitive benchmarking and product failure diagnostics. For these reasons, quantitative measures are now required at virtually every level of manufacturing throughout the automotive industry. Process capability indices were introduced to serve as the language for these quantifying measures (Kane 1986). They have since been widely accepted and adopted in industry and are now often the mandated units of measure for describing process quality. Unfortunately, the current capability indices in use are often inconsistent measures with little 1

or no physical meaning. This is especially evident when one tries to extend the common onedimensional indices to measure the quality of a process with multi-dimensional characteristics and tolerance zones. In this paper we describe some of the current methods used to measure multidimensional process capability, demonstrate how these methods are inconsistent and misleading, and define solution techniques which should be valid for any multi-dimensional characteristics and complex tolerance regions. Our initial motivation was prompted by a question of how to measure capability on the production of transmission gear carriers. A gear carrier is an essential part of an automatic transmission (Figure 1). The primary function of a gear carrier is to hold a number of small pinion gears in a specific location and angular orientation with respect to the centerline of the transmission. Shifting from gear to gear is then accomplished by rotating these small pinion gears, either in conjunction with or independent of other gearing systems within the transmission. Since the entire power of the automobile flows through this gear system, the gear carrier must be held to very tight tolerances for location and angular orientation of the pinion gears. Both of these characteristics are obviously controlled by the location of the holes drilled in the upper and lower tier of the carrier. These upper and lower holes form a coaxial hole pair which determines the eventual location and angular orientation of the shaft for one of the four pinion gears. From a design standpoint, there are three separate two-dimensional tolerance zones for each coaxial hole pair. These three tolerances control the location of the hole on the top tier, the location of the hole on the bottom tier, and the angular orientation of the gear by constraining the location of the bottom hole with respect to the location of the top hole. Many different manufacturing processes are used to produce the coaxial holes in the gear carrier. They vary widely by method and approach but in general they are all responsible for producing the proper locations and orientations of all four coaxial holes. Since quantitative benchmarking and comparison data is needed among the processes, it is imperative that reliable methods exist for quantifying the ability of each process to conform to the tolerance requirements for the four coaxial holes. While there are a variety of reasonable methods for this quantification, we concentrate on those that complement and can be directly related to existing process capability indices such as Cp and Cpk (Kane 1986). Modern manufacturers have invested a great deal of time and money in training 2

engineers, managers, and technicians in the use and interpretation of these indices. In fact, process capability indices have become a familiar, preferred, and often required language for quantifying process performance. For example, Ford's Q-101 manual (1990) states the requirement that all manufacturing processes must achieve minimum process capability indices Cp and Cpk of 1.33 or better (1.0 for previously tooled parts). For the gear carrier process it is common to apply Cp and Cpk to measure-the process capability for a single hole location. However, as illustrated by the examples in the next section, computing Cp and Cpk in the conventional manner for a hole location can easily produce misleading results. Indeed we show that current process capability indices are insufficient to properly quantify the capability for producing even a single hole. Thus they are inappropriate to accurately quantify the capability of the production and tolerance system for the gear carrier. This situation motivated us to develop a new system of capability measurement which can consistently and accurately quantify process quality for a wide variety of multivariate process situations. 2. Current Methods for Bivariate Hole Location Problem A nonconforming part is defined in ANSI (1982) to be: any part that falls outside the tolerance zones. It is therefore clear that the most basic measure of manufacturing process capability is the proportion or fraction of nonconforming parts produced. In this section, we present definitions of the commonly used indices Cp and Cpk, and show that Cpk can be an inconsistent measure of the proportion of nonconforming parts. We explain how Cp and Cpk are currently computed for the bivariate hole'location problem and why they are misleading measures of process performance. These examples motivate the definition of new indices that are applicable to both univariate and multivariate processes. For a univariate process, let T = the design target or nominal, A = the process mean, a= the process standard deviation, U and L = the upper and lower design specifications. 3

In this paper, we assume that the target value T lies in the center of the specification limits. The univariate Cp index is defined in Kane (1986) as U-L CP = (2.1) and is a function of only the process variation. Thus, Cp measures the potential performance of a process: that which could be attained if the process mean is equal to T. Note that for a given value of Cp, the potential proportion of nonconforming parts can be directly determined. In particular, the associated proportion nonconforming (when the univariate process is normally distributed) is p = 2[1 - (3Cp)] (2.2) where $(x) is the cumulative standardized normal distribution. The Cpk index, related to the Cp index, is a measure of the actual performance of a process with the process mean not necessarily equal to T. Again, from Kane (1986), Cpk = Cp(l - k) where k- 2 T- (2.3) U-L Note that for a given value of Cp and k, the actual proportion nonconforming is given by p = 2 -'(3Cp(1- k)) - I(3Cp(1 + k)). -(2.4) The index k, which compares the difference between the process mean and the target to half the allowable spread, is widely used in Japanese industry as a measure of process centering (Kane 1986). Both Cp and Cpk are generally used to quantify process performance under the conditions of normality or near-normality of the realized part dimension (Pearn et al. 1992). 2.1 Application of Cp and Cpk to the Hole Location Problem In order to demonstrate their deficiencies, we now present (via examples) commonly used methods that reduce the inherently two-dimensional problem to one-dimension so that the traditional capability indices can be computed. Consider a manufacturing process which is responsible for drilling 4

a hole in a plane. A target location (xo, yo) for the center of the hole is specified where x0 and yo are measured from some fixed reference point on the plane. The specification region for these characteristics is frequently in the form of a circular or elliptical region (ANSI 1982). It follows that the hole location process relative to the fixed reference point on the plane is a two-dimensional random variable (X, Y). We assume that (X, Y) has a bivariate normal distribution with process mean t = (pA1,y) and variance-covariance matrix ~. Furthermore, we assume the specification region is a circle with a center of (zo, yo) and radius r. Note that in this paper when we refer to "hole location", we mean the location of the center of the produced hole. It is assumed that the radius of the produced hole is separately toleranced. Example 1 The true position deviation for the location of the center of the hole is defined to be: D= /(X-o) + (Y-yo)2, where (xo, yo) is the target location, and (X, Y) is the realized hole location, both relative to some fixed reference point on the plane. Using D as a one-dimensional proxy for the two-dimensions of (X,Y), Cp and Cpk can be computed by using the mean and standard deviation of D. One immediate weakness of this method is that D, being strictly nonnegative, is not likely to be modeled well by a normal distribution. Hence the usual interpretation of Cp, and the well known statistical methods to estimate Cp based on independent and identically distributed normal samples no longer apply. Consider a process A that produces a hole location process (XA, YA) that is bivariate normal with ILA = (o0, y0) and 2- 2 rc S~~A = r2c2 where c is a positive constant. We shall compare process A to another process B that produces a hole location (XB, YB) which 5

is bivariate normal with JUB = (xo - r/2, yo) and 0 62 where 6 is a positive constant. The contour plots for the respective density functions are given in Figure 2 for c = 0.75 and 6 = r/36. Let DA = /(XA - o) + (YA yo) and DB = V/(XB - )2 + (YB o)2 - It can be verified that VAR(DA)= 2 ) VAR(DB) VAR(XB) = (r/6)2, when 6 << r/6. Thus VAR(DA) < VAR(DB) if c < 1/(2V2 - 7r/2) P 0.76 and so if 6 and c are small enough, the Cp of process A is greater than the Cp of process B. This result suggests that the potential capability of process A is greater than the potential capability of process B. However, a visual examination of Figure 2 indicates that process B has a better potential for producing fewer nonconforming parts if the process mean can be centered at the target. Example 2 The part center radial deviation measures the radial deviation from the realized point to a reference point on the plane outside of the specification region. Generally, this reference point is the part center or a locating pin position on the part. If this reference point is (xo, YO) then the part center radial deviation is given by D'- = / (X - x~) + (Y - y)2. Since the expected value of D' is usually much larger than its variance, the probability distribution of D' can be usually be approximated fairly well by a univariate normal distribution. 6

Let us compare process E that produces (XE,YE) and process F that produces (XF,YF). Assume 4o = 0, yA = 0 and let R be the distance between the reference point (0,0) and the center of the specification region of radius r. The random vector (XE, YE) is assumed to follow a bivariate normal with AE = (0, R) and - 0o] SE = 1 2' 0 a0r and (X, YF) is bivariate normal with /F = (0, R) and a22 0 F =. 0 al where al >> a2. An example of this situation is illustrated in Figure 3. As can be seen in this figure, the distributions of the two processes are identical except for a 90 degree rotation of the principal axis. Regarding the production of nonconforming parts, process E and process F are equivalent. However, the part center radial deviation method will report that process E has a much greater potential capability than process F. This is because DE= + YS and DF = v/X + Y2. and the computed values of C, for processes E and F are proportional to the standard deviation of D' and Dp respectively. It can be seen that the standard deviation of D' is less than that of D'by the following argument. If al is held constant, and Z is a standard normal variate, E(D 2) = E(DF2) = a2 + a2 + R2 and as a2 goes to zero E(DE) = E[ rZ2 +R2] but E(D) = R as a2 goes to zero. Hence provided that a2 is small enough, it is clear that VAR(DI') < VAR(DI). We have shown that computing process capability for two-dimensional characteristics using either of these methods can produce inconsistent and inappropriate results. In the next example, we show that Cpk can be a misleading measure even in one dimension. 7

Example 3 Consider two processes G and H where process G has Cp = 1.01 and Cpk = 1.01 and process H has Cp = 1.33 and Cpk = 0.97. Then using a criterion such as that stated in Q-101 (Cpk must be > 1.0), process H would be stopped and 100% inspected while process G would be left to run. However by Equation (2.4), process G produces an actual proportion of p = 2 - 4(3.03) - P(3.03) = 2.446 x 10-3 nonconforming parts, while process H produces an actual proportion of p = 2 - $(5.07) - 4(2.91) = 1.807 x 10-3 nonconforming parts. Therefore, basing a control decision on Cpk can lead to misleading and potentially costly results. This example simply highlights the fact, stressed by several authors (Boyles 1991), that Cpk is not a consistent measure of the actual proportion nonconforming. Specifically, given any value of Cpk, the actual proportion nonconforming is bounded above and below by 2[1 - (3Cpk)] and [1- $(3Cpk)]. Since Cpk behaves inconsistently in one-dimension, it is not surprising that it does even worse in two-dimensions. 3. General Criteria for Evaluating Process Performance There has been a considerable amount of disagreement and debate regarding the criteria that should be used in quantifying the capability of a production process to meet its design tolerances. Much of this debate has focused on comparing the merits of various individual indices. However, it is clear that even -for a simple one-dimensional process characteristic no one index can adequately describe all relevant process information. Other authors have introduced capability indices for multivariate processes. For example, Chan et al. (1991) introduced the Cp,, index and Pearn et al. (1992) defined a multivariate generalization of Cp and Cp,, which they called the,C2 and vC2 indices. However, since each of these indices maps quality into an index value differently, it is difficult to make consistent interpretations when using them. We believe more effort needs to be devoted to tying together various measures to form a complementary and comprehensive set of indices to describe all useful process information. A capability report should accurately and consistently summarize process information to assist engineers and managers in decision making. Our objectives are therefore to report measures that have direct physical meaning and supply information relevant to process description and improvement. Specifically, we recommend that three separate criteria should be used to judge or evaluate the performance of a manufacturing process. 8

1. p = Actual Proportion Nonconforming The actual proportion nonconforming currently being produced by the process, p, is a function of both the current process mean and process variation. It reflects the current state of the process and assumes that no actions are taken to shift the process mean or to reduce the process variation. 2. p* = Potential Proportion Nonconforming p* measures the minimum possible proportion of nonconforming parts that can be achieved through a simple shift of the process mean. In the common case of a symmetric distribution (such as a normal), this minimum value will be attained when the process mean is centered within the tolerance zone. Since most short term engineering procedures for optimizing process performance are limited to a shift of the process mean, this measure conveniently summarizes a best case scenario for potential proportion nonconforming without requiring a reduction in the process variation. Strictly a function of the process variation, the potential proportion nonconforming can be decreased only through a reduction in process variation. 3. k =Process Centering Index For the symmetric case, when the actual and potential proportion nonconforming are not equal, the process is not perfectly centered. However, it is generally unclear precisely how poorly centered the mean is with respect to target. Therefore a natural quantitative measure is to compute the ratio of the difference between the process mean and the design target to the allowable tolerance. The k index defined in Equation (2.3) is an example of such a measure for a one-dimensional process. Note that this measure does not depend on the process variance. Besides the three criteria presented above, a fourth possible criterion that may be useful in describing process performance is one that incorporates a general loss function approach (e.g., see Alkhairy and Staelin 1992, Phadke 1988, and Taguchi 1989). In the case of a 0-1 loss function, the actual proportion nonconforming is equivalent to the expected loss. Many authors have convincingly argued for the use of an alternative loss function, such as one that is quadratic in the deviation from target values. Although such a measure might be of eventual interest, our concern here is with developing indices that are a) direct extensions of existing univariate indices to multivariate situations; and b) easily interpreted and adopted in practice. 9

3.1 Definitions of Cp* and Cp using p and p* While it is informative to report p and p*, they are not directly compatible with most current corporate requirements. It is therefore convenient to transform p and p* into equivalent indices that retain their physical meaning (actual and potential proportion nonconforming) but allow communication in a language familiar to engineers, managers and technicians. In particular, the transformation we use ensures that when the underlying process follows a univariate normal distribution, the new capability index Cp* to measure the potential proportion nonconforming is the same as the Cp index defined in Equation (2.1). This is of practical importance since Cp has become a well known, understood, and deeply entrenched fixture in industry for measuring potential capability. Consider the function f(xz)= ~-1 (1- ). (3.1) A plot of f(x) is given in Figure 4 where x is in logarithmic scale. Define a new index of the actual proportion nonconforming, Cpp, to be Cpp - f(p). (3.2) Note that Equation (3.2) is the equivalent to p = 2[1-,(3pp)]. The actual capability of a process is then obtained in the following sequence: First a value of p, the actual proportion nonconforming, is somehow obtained (this can be done analytically, by simulation, statistically, etc.). Then, using Equation (3.2), the value of Cpp is reported out for compliance with corporate requirements. Similarly, we define Cp* by means of the potential proportion nonconforming p*: Cp. - f(p*). (3.3) where, again, first p* is obtained (measured, estimated, calculated) and then Cp* is reported out. We see that Equation (3.3) can be written p* = 2[1- -(3Cp*)], which is identical to the traditional definition given by Equation (2.2). 10

3.2 Definition of kL a Process Centering Index for Bivariate Location To illustrate the use of Cp* and Cpp as defined in Section 3, and to introduce a generalization of the centering index k, consider a manufacturing process which drills a single hole in a plane. The tolerance region for the hole is defined as the set of all (x, y) pairs such that (x-x)2 + (- 1.'(3.4) a2 2 This (elliptical) region is shown in Figure 5 for the case where (xo, yo) = (0,0 ). (Such an elliptical tolerance zone for a hole center is appropriate where deviation in one direction affects product functionality more than deviation in the other.) A bivariate equivalent of the process centering index. k can be defined by noting that the k index defined by Equation (2.3) is the ratio of the distance that the process mean is off-center with respect to the length of the tolerance zone. Thus, for a bivariate characteristic and an elliptical tolerance zone centered at (xo, yo), we generalize the concept of "distance" so that ( lf^A -o ) 2 I( ^y- - YO 2 kL + b (3.5) where (,ax, ^,) is the process mean of the bivariate characteristic. Note that both k and kL reflect only the amount that the process mean is off-center. They provide no information about the direction of this deviation with respect to the target. 4. Capability Analysis for a Single Hole Location Problem 4.1 Computing Actual and Potential Proportion Nonconforming In this subsection, we present an efficient method for computing p and p* for the specific case of a bivariate normal characteristic and an elliptical tolerance zone. In particular, let (X', Y') represent the hole location characteristic of the process. Assume that (X',Y') follows a bivariate normal distribution (BVN(t', E')) with mean tL = (4/' R) and variance-covariance matrix A', where 11 P. I 11

(In practice, iL' and E' are of course unknown and must be estimated. If we use the sample mean A' and sample variance-covariance matrix E, then the method described below will produce the corresponding estimates for p and p*. Computing the accuracy of these estimates, and associated error bounds, is discussed later.) The probability that (X',Y') will fall inside the elliptical region of Equation (3.4) can be computed by first transforming the coordinate system so that the region is circular with radius r. To do this, we transform the BVN(1', Z') to a BVN(i- = (C,Ay), S) where Ua2 ~ X = O L o j [o~"] It can be shown that r =b, b/, cos0, b2a2 2 COS2 6 2bpa'a' cos 0 sin 2. aa= -+ a- + aT- sin2, b2a sin2 0 2bpa'a cos 0 sin 0 2 c2 a = - + a/ cos 2, a2 a where 2pbal a' Y a/1^ a ya a 2 - tan0 = 2pboc{U z < ( < () The above transformation ensures that (a, p'y) lies in the first quadrant of the new coordinate system, and that a. > ay. Furthermore, the elliptical region of Equation (3.4) becomes the circular region x2 + y2 < r2 Using a formula from Groenewoud et al. (1967), the probability that a bivariate normal random variable falls outside a circular region of radius r is given by P=l-aZESm (pff2) > Dm (X'-+ ^ ) (4.1) Jsn2 2, 2 2 4 2 4 where 2 2 ] is the greatest integer less than or equal to 2 X7 a is given by a { ( 2 12

and Smn(x) is the Poisson tail probability m-1 k Sm() = 1- e- Z' k=O The terms D (u, w, y) can be computed recursively from D~(u,,y) = 1, D~(u,w, y) = y+ u Do(u, w,y) = 1 m- u +y] D~(u,w,y)-uyD (u,, ) for m >3, (4.2) D( M- w y iU WIjUl.) D(u, w,y) = uD _(u, w,y)+ _Dl_2(u, w,y) for j 1 and m > 2j + 1. (4.3) The initial values for D2(u, w, y) and D~(u, w, y), along with Equations (4.2) and (4.3), completely define the Di term for all relevant ranges of the summation in Equation (4.1). For computational purposes it is necessary to truncate the summation over m at some value m = N. Groenewoud et al. (1967) show that when this is done, the error in the computed probability, EN, has the bound p 32N(N + 1) r EN "< e-. 1 - + 2) where P 2 c - To ensure that EN < K, it is possible to find the smallest value of N such that e-2 /32N(N + 1) < K (4.4) N!(N+1-1P2)- This integration method provides computations for any practical error EN in a matter of seconds on simple desktop PCs. Note that the above algorithm is used twice: to determine both a) the actual proportion nonconforming p, and b) the potential proportion nonconforming p* (by setting the process mean at the design target, I = ly = 0). 4.2 Computing Cp* and Cpp Indices Once p and p* have been determined, we can compute Cp* = f(p*) and Cpp= f(p). This requires a computational method for accurately evaluating the tail area of the standardized normal curve V L e-z2/2dx = 1 - (Z). 13

Hill and Joyce (1976) present a fast and accurate computational algorithm for this computation making use of two significant results from Abramowitz and Stegun (1972). If z lies in the central area of the normal curve, then the algorithmic method used is'the convergent series 2 /2 Z X2 /2, Z3 Z5 Z7 f2 I -a;/2 e d = z- + -- - + Jez2 e-32/2dz+ —- 3 x5 3x5x7 If z lies in one of the tails of the normal curve, then the method used is the continued fraction ez2/2 e-2/2 dx = 1 z+ 2 z+ z+ They recommend the changeover point between the two methods to be lzJ = 2.32. This procedure is extremely fast and computes tail probabilities to virtually the accuracy of the machine for z < 7 and to one digit less than the accuracy of the machine for z > 7. In order to find Cp*, it is necessary to solve for Cp* in * = 2[1- 4(Cp*)] = 2((-Cp*). (4.5) A bisection root finding method can be used to solve for any p* calculated from Equation (4.1). To solve for Cpp, p* is replaced by p in Equation (4.5). In practice, implementation of this method is limited only by the length of the double precision variables of the computer. Although this is both machine and compiler dependent, 19 digits of precision for floating point variables on a desktop PC is not uncommon. Therefore, numerical errors only begin to occur for Cp* (or Cpp) values greater than or equal to approximately 3.0, a value that is essentially indistinguishable from higher index values since it implies that the proportion nonconforming is on the order of 10-19. 4.3 Example Consider the scatterplot presented with Table 1, plotted over a circular tolerance zone with radius r = 0.1 mm and target value of (0,44.45) mm. (The data is from an actual manufacturing process: we do not include the raw data set for proprietary reasons). Using the original data we have 14

t' = (0.0042,44.4667) and =' [ 5.83 2.47 x 1-4. 2.47 2.58 Applying the coordinate transformation gives 6 = 28.32 degrees, ft = (0.0116,0.0127) and 7.15 0 E = X 10-4. 0 1.38 To achieve a desired accuracy of K = 1.0 x 10-12 in computation of p and p*, Equation (4.4) produces a cut-off value of N = 87. Using Equation (4.1) then gives p* = 2.296 x 10-4 and p = 6.689 x 10-4 from which Equation (4.5) gives Cp* = 1.23 and Cp = 1.13. Using the sample mean estimate in Equation (3.5) gives kL = 0.17. In order to make statistical confidence statements about the true process values of Cp*, Cpp, and k from the sample data, a bootstrap simulation was performed. The bootstrap method resamples (with replacement) from the original data set and computes estimates of Cp*, Cpp, and k for each run. The result is an empirical distribution for each of these sample estimates. From these distributions, confidence intervals on the true unknown parameters are generated. See Efron (1982) for more details on the bootstrap and other resampling procedures. For this example, 10,000 bootstrap replications were performed. The descriptive statistics for the resulting Cp*, Cpp, and k are presented in Table 1. From this table, we see that the 95% confidence interval bounds are: 1.11 <Cp* < 1.41 (23.4ppm < p* < 868ppm), 1.00 < Cpp < 1.33 (66.1ppm < p < 2700ppm), 0.13 < kL < 0.22. The general usefulness of kL is demonstrated in this example. Suppose the minimum possible location adjustment for this process (due to limitations on machine control) is 10% of the tolerance radius. In this situation therefore, any kL < 0.10 indicates that the centering of the process cannot be improved. Since the 95% lower bound on kL is 0.13, we can conclude that the process is not centered and could be improved through a mean shift. The estimate of kL = 0.17 immediately indicates that the process mean is estimated to be off-target by 17% of the allowable tolerance. 15

Note that kL does not by itself provide directional information and therefore care must be taken in interpreting hypothesis tests. As shown in Table 1, in this example of the 10,000 bootstrap resamples taken, 9408 of the sample means were in the positive-x, positive-y quadrant. The remaining 592 sample means were in the negative-x, positive-y quadrant. This additional information supports an inference that the mean of y is biased in the positive direction, and the mean of x is also biased in the positive direction but to a smaller degree. 5. Capability Analysis for a Coaxial Hole While single hole location tolerances are more commonly encountered in design specifications, coaxial hole tolerances are also used for critical dimensions. We now show how our method can be used for measuring process capability for a coaxial hole tolerance system. A two-tier coaxial hole usually constrains a pin or shaft to lie in a proper location and angular orientation. For a single hole, there is one two-dimensional tolerance region for the location of the hole center. Since a two-tier coaxial hole has two separate holes, there is a pair of two-dimensional tolerance regions, one each for the location of the top and bottom holes. As discussed in Section 4, these tolerance regions are usually circular or elliptical. Hereafter, we refer to both of these tolerance regions as location tolerances. A further constraint is often necessary to control the angular orientation of the shaft. This is usually expressed as a two-dimensional tolerance region for the relative location of the coordinates of the center of the bottom hole with respect to those of the top hole. (By convention, we define the top hole to be the first hole produced by the tool.) This is justified since in most common process plans, once the tool has produced the top hole, it proceeds down to produce the bottom hole. In this case, the location of the center of the bottom hole with respect to that of the top hole reflects the drift of the tool after machining has begun. Hereafter, we refer to the tolerance which constrains this drift as the angular tolerance. While the two location tolerances are fixed with respect to the target location of the hole centers, the angular tolerance depends on the realized location of the top hole for a given part. For example, if (xT, YT) and (XB,YB) represent the realized location of the centers of the top and bottom holes respectively and the specified angular tolerance is circular, then the angular tolerance is a circular 16

region centered at (xT, YT)The overall tolerance region for a coaxial hole is therefore defined to be the intersection of these three two-dimensional tolerances (see Figure 6). Any produced two-tier coaxial hole which fails to meet all three tolerances simultaneously is considered, by ANSI definition, to be nonconforming. 5.1 Definition of Capability Measures for Coaxial Holes Computing p and p*, the actual and potential proportion nonconforming for the coaxial hole tolerance system, shown in Figure 6, requires finding the probability that a part fails to simultaneously conform to all three tolerance requirements. An appropriate measure of process centering is slightly less straightforward. In particular, it seems pointless to have a single scalar index of a process that includes two separate notions of centering: a) the location of a mean circle center relative to a fixed reference point; and b) the location of a second mean circle center relative to the produced first center. However, the process centering capability for a coaxial hole could be measured by two separate k indices: kL, the ratio of the deviation of the top hole process mean from its target; and kA, the a ratio of the distance between the process means of the top and bottom hole with respect to the width of the angular tolerance region. These can be formally defined as follows: Let (rxi, Lyl) and (,zx2, sy2) be the process mean of the top and bottom hole location respectively. If the location tolerance regions and angular tolerance regions are both circular, with radii rL and TA respectively, then. we define k -V/(/1 - x0)2 + (Ati - yO )2 (5.1) kL = (5.1) and k_ I(1iL - r2) + (11y - 2)2 (5.2) kA = (5.2) TA where (xo, yo) is the target location of the top hole. Note that kL reflects the ability of the manufacturing process to accurately position the drilling tool at the proper location prior to machining. The index kA reflects how well the drill maintains its angular orientation perpendicular to the drilling surface throughout the machining process. Since both of these indices measure independent physical properties, reporting both of them is consistent with the overall goal of summarizing all relevant process information in an easily understood language. 17

5.2 Computation of Capability Indices for Coaxial Holes Recall that (XTYT) and (XB, YB) represent the (random variable) locations of the top and bottom holes. Assume that X = (XT,YT,XB,YB) follows a four-dimensional normal distribution with mean IZc = (Aili,!Aly, x2z, lz2y) and variance-covariance matrix Ec. Computation of p and p* requires a method for integrating this four-dimensional normal distribution over the intersection of three distinct tolerance regions. Monte Carlo integration is often useful, although time-consuming, for computing p and p* since it can integrate high-dimensional distributions over any complex region. In the example below, Monte Carlo integration, programmed on a Risc 6000 IBM workstation, generates approximately one million four-dimensional random variates every 3 minutes. An accurate estimate of a p (or p*) for a proportion nonconforming below 1 part per million, needing at least 50 simulated non-conforming observations for statistical significance, may require an hour or more. We are currently researching methods to carry out this integration efficiently. Once p and p* have been computed, Cpp and Cp* can be obtained by solving Equation (4.5) numerically, as discussed in Section 4. 5.3 Example This method is best demonstrated through the use of an example collected from actual manufacturing data. A gear carrier has four different coaxial hole pairs with top and bottom hole location target values of (44.45,0), (0,44.45), (-44.45,0), and (0,-44.45) mm. A circular tolerance region with radius of 0.1 mm specifies the location of the top and bottom holes. The angular tolerance is specified as a circular region with a radius of 0.075 mm centered at the realized hole center of the top hole. (Angularity is usually specified with a tighter tolerance than location since the angle of the pinion shaft is a more critical element for proper gear function). 78 sample observations were collected, each consisting of eight hole (center) location coordinates. This data is shown in Figure 7. The top hole for coaxial hole. 1 was used as a location reference by the coordinate measuring machine, so the sample data for this y-coordinate is always zero. Figure 7 also shows the angular tolerance regions, as well as points obtained by subtracting the bottom hole coordinates from the top hole coordinates, since the target value for each of the four angular tolerances is (0,0). 18

5.4 Capability for Each Coaxial Hole Pair We first concern ourself with the computation of p and p* for each of the four coaxial hole pairs. For each coaxial hole pair, the capability of the top and bottom hole can be computed independently using the method given in Section 4. The capability of the angular tolerance and hence the capability of the entire system of tolerances must be computed using the coaxial hole method described above. To evaluate the four-dimensional normal integrals as discussed in Section 5.2, Ac, and Ec are estimated for each coaxial hole pair, and Cholesky decomposition is applied to find a lower triangular matrix A such that ec = AAt. The resulting lic, Sc and A for each of the four coaxial holes are shown in Table 2. Next, four independent N(0,1) random variates Z1, Z2, Z3 and Z4 are generated, and X = Ac+ZAt is computed, where Z = (Z1, Z2, Z3, Z4). The random variate X is thus normal with mean AiC and variance-covariance matrix c. Each generated variate X is checked to see if it satisfies the tolerance system shown in Figure 6, i.e. "that it conforms". The proportion nonconforming is determined by dividing the observed number of nonconforming items by the total number of fourdimensional random variat'es generated. This same random variate can also be used to estimate the potential proportion nonconforming by centering the process mean vector at the design target. The results, along with Cp*, Cpp, and the kL indices for both top and bottom holes and the kA index are shown in Table 3. Note that these indices for the coaxial hole system complement, but do not. replace, the information provided by the indices for each single hole. As in the single hole case, a bootstrap analysis was performed on the coaxial data. The resulting 95% confidence intervals for the individual tolerances and the entire coaxial hole are presented in Table 4. These results are based on runs of 10,000 for the top and bottom hole and 1000 runs for the angular and entire coaxial holes. 6. Capability Analysis for the Entire Gear Carrier The analysis in Section 5 produced capability measurements for each individual coaxial hole. An interesting and challenging issue is the determination of the capability of the system comprising all 19

four coaxial holes pairs. Since a single process is responsible for making all four coaxial holes, the performance of the process over the entire tolerance system is relevant. The major difficulty with computing the capability of the entire system is that it is a 16dimensional process (eight holes, each with two dimensions describing its center), involving a total of 152 unknown parameters that must be estimated (16 components of the mean vector and 136 components of the symmetric variance-covariance matrix). This is clearly infeasible unless a large data set is available. Furthermore, it would be difficult to test for assumptions about statistical distributions. Although these challenges are substantial, it is likely that a large number of the elements in the 16-dimensional variance-covariance matrix are approximately zero. (We have collected about a dozen process samples and the empirical evidence supports this position.) Assuming normality, such zero correlation allows a decomposition of the 16 variables into smaller sets of independent variables, whose individual conformance probabilities can be multiplied. However, even without making these assumptions of independence, there are simple bounds on p for the entire four hole coaxial system that can be computed. Let pi represent the actual proportion nonconforming for each coaxial hole i, i = 1,2,3,4. A gear carrier is defined to be nonconforming whenever any of the coaxial holes fails to simultaneously meet the three tolerances as described in Section 5. Then it can be shown that bounds on p for the entire system are given by PL max(pl, p2, p3,4) < p < min(l, pi + P2 + P3 + 4) PU. (5.3) It follows that f(pu) < Cpp < f(PL), where f is the transformation defined in Equation (3.1) of Section 3. The same bounding inequalities clearly also hold for p* and Cp*. For the data of -Figure 7, this results in 0.76 < Cp* < 0.89 (7318ppm < p* < 22,966ppm), 0.61< Cpp < 0.74 (26,014ppm < p < 65,475ppm). Process centering for the entire four hole system can also be obtained by generalizing Equations (5.1) and (5.2). Let ( 21 iy1li) and (IL:2i, yt2i) be the process mean for the top and bottom hole of coaxial hole i, i'= 1,2,3,4. As in the coaxial hole case, we define two separate k indices, one for the overall location of the drill pattern and the second for the overall angular orientation. Again, let rL and rA be the radii of the location and angular tolerance regions of an individual coaxial hole system. If (zoc, /Yoc) is the target value for the center of the top hole pattern, then we 20

can define \ 1(4 E lxii - XoC + )L ( i - Yoc i-=1 i=1 kL = rL and kA4 k(_zli -!L2i) + [ E(U/li - 12i)], L i=J LP ~i J= kA = i= i=l rA Estimates kL and kA, for the data in Figure 7 are (kL, kA) = (0.06,0.25). The bootstrap simulation carried out for the coaxial hole also provides sampling information for the entire carrier. For each of the 1000 resampling runs, the bounds of Equation (5.3) were recomputed. From these we can produce "95% bounds " by solving for the interval that is completely covered by 95% of the resampled bounds. From this we form the following 95% confidence interval on the true capability of the entire gear carrier: 0.58'< Cp, < 1.02 (2213ppm < p* < 83,451ppm), 0.47 < Cpp <0.84 (11,537ppm < p < 158,540ppm), 0.02 < kL < 0.10, 0.21 < kA < 0.29. Note that these bounds are fairly wide, and only provide an order of magnitude estimate of the proportion nonconforming. A more accurate analysis would require a larger sample size of the original data. Also note that kL = 0.06 for the entire process while kL (from Table 3) for the individual coaxial holes varies from 0.15 to 0.29. The reason for this can be seen directly from Figure 7 since the average location for the top of each coaxial hole has a greater radial deviation from center than specified. However, these deviations cancel out for the entire four hole system. Thus, the overall k indices have a special relevance when the process drill bits cannot be moved with respect to one another since kL and kA reflect the overall centering (both for location and angularity) of the entire drill pattern. 21

7. Summary and Conclusions The capability indices described in this paper are general and have application beyond the specific gear carrier example with circular location and angular tolerances. The criteria and the transformed indices Cp* and Cpp defined in Section 3 are useful performance measures for any complex multivariate processes. In this paper, we concentrate on circular tolerance regions since hole location tolerances are usually circular. However, the methods developed here have obvious extensions to rectangular or other general shape tolerance regions. Although we only discuss how the actual and potential proportion nonconforming can be computed under a multivariate normal assumption, extension to other statistical distributions is clear. Furthermore the assumption of multivariate normality is not unreasonable since it is likely that, after suitable transformation of one or more of the variates, data can be considered to be approximately normal. Koziol (1986) presents a comprehensive review of available methods for assessing multivariate normality. While the capability indices described in this paper have meaningful physical interpretations, statistical properties of their estimators are not yet fully explored. It is unlikely that the distributions of, p* or k can be explicitly determined for the multivariate cases. While other multivariate capability indices with better defined statistical properties have been proposed, these indices have no obvious physical meaning and therefore have the potential for being misunderstood and misinterpreted by engineers and managers. Practitioners are therefore left to choose between indices having nice statistical properties but inconsistent meaning, and those (as presented here) with direct physical meaning but difficult statistical properties. Considering the availability of computer-driven statistical estimation procedures (such as the bootstrap), the choice becomes clear. Acknowledgements The work reported here was supported in part by the National Science Foundation under grant DMC-8903029. This work was also supported by a grant from Ford Motor Company. References [1] Abramowitz, M., and Stegun, I. A. (1972), Handbook of Mathematical Functions, Dover Edition. 22

[2] Alkhairy, A. and Staelin, D. (1992), "The Choice of Quality Cost Function for Product and Manufacturing Design," MIT Department of Mechanical Engineering. [3] ANSI Standard Y14.5M. Dimensioning and Tolerancing, ASME (1982). [4] Boyles, R. A. (1991), "The Taguchi Capability Index," Journal of Quality Technology, 23, 1, 17-26. [5] Chan, L. K., Cheng, S. W., and Spiring, F. A. (1991), "A Multivariate Measure of Process Capability," Journal of Modeling and Simulation, 11, 1-6. [6] Efron, B. (1982), The Jackknife, the Bootstrap and Other Resampling Plans, CBMS-NSF Regional Conference series in Applied Mathematics. [7] Ford Worldwide Quality System Standard Q-101, Ford Motor Company (1990). [8] Collected Algorithms from ACM, Association for Computing Machinery (1976). [9] Groenewoud, C., Hoaglin, D. C., Vitalis, A. (1967), Bivariate Normal Offset Circle, Probability Tables with Offset Ellipse Transformations, Vol. 1, Cornell Aeronautical Laboratory of Cornell University. [10] Kane, V. (1986), "Process Capability Indices," Journal of Quality Technology, 18, 1, 41-52. [11] Koziol, J. A. (1986), "Assessing Multivariate Normality: A Compendium," Commun. Statist. Theor. Meth., 15, 9, 2763-2783. [12] Pearn, W. L., Kotz, S., and Johnson, N. L. (1992), "Distributional and Inferential Properties of Process Capability Indices," Journal of Quality Technology, 24, 4, 216-231. [13] Phadke, M. S. (1988), Quality Engineering Using Robust Design, Prentice Hall. [14] Taguchi, G. (1989), Quality Engineering in Production Systems, McGraw Hill. 23

Top Hole Located 01- 1 Appropriately Gear Axis Parallel to Part Center Axis Bottom Hole Located Appropriately Figure 1: A Four Hole Gear Carrier with Three Tolerance Requirements per Coaxial Hole

Process APrcs Tolerance ~"'; =~''::~:~:~:~:~: I~~::Z one.........I.f C~5:~..........~ ~ ~ Z~i~~~.......... sf3=:. ~:2::::;;::............ ~ ~ ~ r: ~I~rf::~:i:~ ~cl~;~;;r;:~:~~:2.i~~'::r:~:~~;.~I~~;s~:;.;t~~:::::01: YO) Figure: Contors of he pdf or HoleLocatins Procsses A nd B. he TruePositio Radia Deviaton Methd willReport hat Pr A is Potentially More Capable than Process B.s~'r~.~;~+ztr

Process C Process D R ~~~~~~R I (0,10) Part Center (0,0) Part Center Figure 3: Contours of the pdf for Hole Locations Processes C and D. The Part Center Radial Deviation Method will Report that Process C is Potentially More Capable than Process D.

v, o o -=. *- =o 0 o —-------- o 1m *oE ~ ttmmm tmm mmm 0 4 a- X r - ^ -"'*^~- - - " ^ -- y ^- --------- O - - - Q - -........ i.....:.......- -. - 0 0 "_"" _/ 0;=:::^=- ------ - I = ======= 2 =, - - -, o —S* - 7^c-_^. —^~-^ — -^-^_^_ * ( ----- - -- - - --------- - ---- i 1 _... _ o. —... =L== —-- ======== —--, ~- --- V- ~, t-/ - -= -- - - - _>I,:' ": _ ---------- --.-;c IJ _ _1 _*~, 3 s e on" U FZ~~~~~~~~~~~~~~~~~~~~~~~~T

y Toleranc Zon xFigure 5: Elliptical Tolerance Region in (z, y) Plane b Figure 5: Elliptical Tolerance Region in (x, y) Plane

Elevation View Top Hole Location """"""' Actual Top Hole Tolerance -.:.. Location for an Observed Part.... __ [ Angular Tolerace Bottom Hole - region extended Location - from the Observed Tolerance Top Hole Location Overhead View Top Hole Location Bottom Hole Location Tolerance Tolerance Simultaneous Bottom Hole Tolerance Realized Angular / Tolerance for an Observed Part Figure 6: Simultaneous Tolerance System for a Coaxial Hole Consisting of Three Circular Tolerance Regions

Top Hole Tolerance _. I\. _ Angular Tolerance Coaxial.:. I Hole #2 Bottom Hole Tolerance I I Top Hole Tolerance Top Hole Tolerance |ngul^ To Hole Tolerance ~' IB, ~. Hole #3 Hole #1 t ^~~~ —— I~ ACenter Bottom Hole Tolerance I Bottom Hole Tolerance Top Hole Tolerance dJ J Angular Tolerance Coaxial (-^ —~~ ~Hole#4 I Bottom Hole Tolerance Figure 7: Data from a Gear Carrier Sample of Size 78. Note that the Data Plotted versus the Angular Tolerance is Computed by Subtracting the Location of the Bottom Hole from the Location of the Top Hole for each Part

Hole 2 - Top Cp* Cpp k Estimate 1.228 1.134 0.173 Mean 1.241 1.145 0.174 Standard Deviation 0.078 0.085 0.022 Skewness 0.437 0.441 0.253 y Kurtosis 0.391 0.372 0.121 Minimum 0.995 0.885 0.089 44.50 Maximum 1.649 1.625 0.276 Count 10000 10000 10000 95.0 Percentile 1.379 1.294 0.213 5.0 Percentile 1.124 1.016 0.139 4445 44.45 - 97.5 Percentile 1.406 1.328 0.221 2.5 Percentile 1.106 0.996 0.134 99.5 Percentile 1.473 1.401 0.236 0.5 Percentile 1.066 0.955 0.123 44.40 1.0 Percentile 1.081 0.971 0.127 Location of Bootstrap Means by Quadrant No. in Positive-x, Positive-y 9408 - 0 -- ^ - - -' -- - ----- -0.05 0.00 0.05 X No. in Positive-x, Negative-y 0 No. in Negative-x, Negative-y 0 No. in Negative-x, Positive-y 592 Table 1: Scatterplot and Descriptive Statistics for Indices Computed from a Sample of 78 Hole Locations Observations

Coaxial Hole pc,c A 44.479 3.25e- 4 2.11e - 6 2.79e - 4 1.11e - 5 1.80e - 2 0.000 2.41e -6 2.24e -6 5.63e -6 1.17e- 4 1.55e-3 44.471 8.06e -4 -6.89e-5 1.55e- 2 2.74e -4 2.38e-2 -0.014 1.82e - 4 6.14e- 4 3.59e- 3 -3.34e- 3 1.25e-2 0.004 5.83e -4 2.47e - 4 5.12e- 4 1.80e - 4 2.41e - 2 44.467 2.58e-4 2.43e -4 1.62e -4 1.02e - 2 1.24e - 2 -0.008 1.50e-3 1.04e -4 2.12e - 2 2.08e -3 3.24e-2 44.450 2.69e-4 7.47e - 3 6.87e - 3 -2.14e - 3 1.27e-2 -44.469 3.43e -4 -2.89e -5 2.85e- 4 -1.31e- 5 1.85e - 2 -0.010 1.05e- 3 -3.42e- 5 9.02e- 4 -1.56e- 3 3.24e -2 -44.469 8.02e- 4 -1.04e- 4 1.54e-2 -3.16e-4 2.38e-2 -0.028 9.27e- 4 -7.08e- 4 2.78e- 2 -3.52e - 3 1.19e 0.007 5.98e- 4 -2.64e- 4 5.66e- 4 -1.79e - 4 2.44e - 2 -44.463 2.78e- 4 -2.78e- 4 1.84e- 4 -1.08e - 2 1.27e- 2 0.007 1.12e- 3 -2.69e -4 2.32e- 2 -2.23e- 3 2.40e- 2 -44.485 2.86e - 4 -7.30e - 3 8.30e - 3 -3.36e - 3 1.23- 2 Table 2: Sample Means and Sample Variance-Covariance Matrices for Coaxial Hole Problem

Coaxial Hole Tolerances Cp (f*) Cp (p) index Top 1.84 (0.0381ppm) 1.31 (87.12ppm) 0.29 Bottom 1.15 (566.8ppm) 0.95 (4210ppm) 0.25 1 Angular 1.00 (2780ppm) 0.95 (4210ppm) A = 0.21 Coaxial 0.99 (3037ppm) 0.89 (7791ppm) (k, kA) = (0.29, 0.21) Top 1.23 (229.6ppm) 1.13 (668.9ppm) 0.17 Bottom 1.00 (2727ppm) 0.97 (3478ppm) 0.08 2 Angular 0.98 (3445ppm) 0.88 (8291ppm) kA = 0.28 Coaxial 0.92 (5571ppm) 0.85 (11,054ppm) (k, kA) = (0.17,0.28) Top 1.00 (2727ppm) 0.94 (4893ppm) 0.20 Bottom 0.97 (3649ppm) 0.76 (21,906ppm) 0.34 3 Angular 1.01 (2470ppm) 0.97 (3756ppm) kA = 0.24 Coaxial 0.89 (7318ppm) 0.74 (26,014ppm) (kL, &A) = (0.20,0.24) Top 1.20 (303.8ppm) 1.11 (868.5ppm) 0.15 Bottom 0.95 (4497ppm) 0.81 (15,606ppm) 0.36 4 Angular 0.99 (2864ppm) 0.92 (6051ppm) kA = 0.30 Coaxial 0.90 (6934ppm) 0.77 (20,558ppm) (kL, kA) = (0.15, 0.30) Table 3: Capability Measurements for the Coaxial Holes of Figure 8

Hole Tolerance 2.5 Percentile 97.5 Percentile a 2.5 Percentile 97.5 Percentile a 2.5 Percentile 97.5 Percentile a Top 1.64 2.16 0.132 1.15 1.54 0.099 0.25 0.33 0.021 Bottom 0.96 1.47 0.130 0.79 1.21 0.108 0.20 0.31 0.029 1 Angular 0.80 1.39 0.141 0.77 1.27 0.121 0.18 0.26 0.022 Coaxial 0.79 1.34 0.135 0.72 1.14 0.108 Top 1.11 1.41 0.078 1.00 1.33 0.085 0.13 0.22 0.022 Bottom 0.86 1.22 0.093 0.84 1.17 0.084 0.02 0.15 0.034 2 Angular 0.78 1.29 0.128 0.73 1.11 0.094 0.23 0.33 0.025 Coaxial 0.76 1.12 0.96 0.71 1.00 0.076 Top 0.89 1.22 0.069 0.84 1.07 0.058 0.15 0.27 0.025 Bottom 0.84 1.13 0.075 0.68 0.87 0.050 0.28 0.40 0.032 3 Angular 0.81 1.35 0.130 0.78 1.24 0.119 0.20 0.29 0.020 Coaxial 0.76 1.02 0.069 0.64 0.85 0.050 Top 1.09 1.37 0.071 1.00 1.25 0.063 0.10 0.20 0.027 Bottom 0.82 1.14 0.082 0.69 0.96 0.069 0.32 0.41 0.023 4 Angular 0.80 1.34 0.132 0.73 1.17 0.111 0.26 0.35 0.022 Coaxial 0.76 1.09 0.087 0.65 0.93 0.072 Table 4: Confidence Regions for Capability Measurements for the Coaxial Holes of Figure 8 (obtained by bootstrap analysis: n = 10, 000 for top and bottom hole tolerances, n = 1000 for angular and coaxial tolerances)