Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Studies in applied statistics : ship hull design optimization and endogamous group comparisons Thorne, Anona Elaine 1993

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


831-ubc_1993_fall_thorne_anona.pdf [ 3.32MB ]
JSON: 831-1.0086311.json
JSON-LD: 831-1.0086311-ld.json
RDF/XML (Pretty): 831-1.0086311-rdf.xml
RDF/JSON: 831-1.0086311-rdf.json
Turtle: 831-1.0086311-turtle.txt
N-Triples: 831-1.0086311-rdf-ntriples.txt
Original Record: 831-1.0086311-source.json
Full Text

Full Text

STUDIES IN APPLIED STATISTICS:SHIP HULL DESIGN OPTIMIZATIONANDENDOGAMOUS GROUP COMPARISONSbyANONA ELAINE THORNEB.A., The University of British Columbia, 1991A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF SCIENCEinTHE FACULTY OF GRADUATE STUDIES(Department of Statistics)We accept this thesis as conformingto the required standardTHE UNIVERSITY OF BRITISH COLUMBIAAugust 1993© Anona Elaine Thorne, 1993In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.(SignatureDepartment of ^StatisticsThe University of British ColumbiaVancouver, CanadaDate Auaust 31. 1993DE-6 (2/88)AbstractThis paper contains two case studies in applied statistics.The first study originated from a series of seakeeping experiments on twenty-seven dif-ferent ship hull designs. The experimenters sought a method of predicting performancefor hull designs different from those included in the study. The results of these experi-ments had already been analyzed, but we wished to consider the Taguchi method as analternative method of analysis.In this study, we describe the initial design problem and analyses and provide a criticalevaluation of the Taguchi method; we assess its merits to determine its suitability forthe ship hull design problem. We find the Taguchi method unsatisfactory in general,and inappropriate for the ship hull design problem.The second study derives from a case which came to the Statistical Consulting andResearch Lab (SCARL) in the Department of Statistics at The University of BritishColumbia. Having lost all her original data, our client wished to know if it were possibleto "reconstruct" her data or, in any event, to perform a statistical analysis, using thetables of summary statistics in her possession.She wished to model several different pulmonary functions using various demographicand anthropometric measurements as explanatory variables. For each of the differentmodel variants under consideration, she wanted to see if the same model parameterswould be valid for all endogamous groups, or if some groups were sufficiently differentthat the parameter values would be different for these groups.We were able to develop a new general method for application to her problem; thederivation and application of this methodology are presented in the second study.iiContentsAbstract^ iiTable of Contents^ iiiList of Tables^ viList of Figures^ viiAcknowledgements^ viiiI Ship Hull Design Optimization^ 1Summary^ 11 Introduction^ 22 Initial Analyses^ 33 The Taguchi Method of Parameter Design^ 83.1 Quality and the Loss Function ^  83.2 The Taguchi Approach ^  103.3 Experimental Design  143.4 Analysis ^  174 The Taguchi Method and the Ship Hull Design Problem^20iii5 Criticisms of the Taguchi Method^ 225.1 Experimental Design ^  235.1.1 Economy  235.1.2 Interaction Effects ^  235.1.3 Sequential Experimentation ^  245.2 Analysis ^  255.2.1 Marginal Average Plots ^  255.2.2 Response Function  275.3 General Approach^  306 Conclusions^ 31References for Part I^ 32II Endogamous Group Comparisons^ 34Summary^ 347 Introduction^ 358 Background^ 368.1 The Estimation of Physical Fitness ^  368.2 The Data Analysis Problem ^  37iv9 Sufficient Statistics^ 419.1 Notation  429.2 Calculation of Sufficient Statistics ^  4310 Derivation of the Method of Analysis^ 4710.1 Notation ^  4710.2 Derivation of Estimators ^  4810.2.1 A Restricted Model  5410.3 Calculation of Estimates ^  5410.4 Testing Methodology  5611 Results^ 6012 Conclusions^ 67References for Part II^ 68Appendix I: Splus Routines for Part II^ 69Appendix II: Test Results for Part II^ 72List of Tables1^L9 array for elastomeric connector experiment (BT87, p. 22) ^ 152^L8 array for elastomeric connector experiment (BT87, p. 22) ^ 163^Subject Distribution among Caste Groups ^  374^Measurements Collected ^  385^Pairwise Comparisons for Testing Equality of   576^Variables of Interest ^  617^Model Variants Tested  628^Test Results from Bonferroni Comparisons for Men and Women ^ 639^Test Results from Individual and Bonferroni Comparisons for Men . . . ^ 6510 Test Results from Individual and Bonferroni Comparisons for Women . ^ 6611 p-values from Pairwise Comparisons for Testing Equality of p i 's for Men 7212 p-values from Pairwise Comparisons for Testing Equality of p i 's for Women 7313 Test Results from Bonferroni Method ^  74viList of Figures1 Traditional loss function (Ross, p. 42) ^ 92 Parabolic loss function (Ross, p. 42) 103 Result of Heat Treatment Experiment (Ross, p. 47) ^ 134 Factor levels for elastomeric connector experiment (BT87, p. 21) ^ 145 Crossed array for elastomeric connector experiment (BT87, p. 22)^. . . . 166 Average S/N ratio for each control factor level (BT87, p. 23) ^ 187 Mean responses for each control factor level (BT87, p. 23) ^ 198 Interaction and marginal average plots for the 2 2 design (Ryan, pp. 34-35) 269 Marginal average plots for the revised 2 2 data (Ryan, p. 35) ^ 2610 A simple example with two controllable factors (Hendrix, p. 76) ^ 2811 Results for Variants with Significant Differences ^ 64viiAcknowledgementsI would like to thank James Zidek for his guidance and help throughout the developmentof this thesis and for being an important role model for my approach to the practiceof statistics in general. I would also like to thank Malcolm Greig for suggesting theinnovative solution to the group comparison problem and for providing advice duringthe implementation of the solution and preparation of the manuscript. I am grateful aswell to Michael Schulzer for his careful reading of the manuscript and valuable comments.In addition, I must acknowledge the work of Lakshmi Reddy, who brought the groupcomparison problem to our attention and with whom I had many instructive conversa-tions. I would also like to express my appreciation to Mohan Delampady for his adviceand support during my initial study of the ship hull design problem.Finally, I would like to thank the Natural Science and Engineering Research Council ofCanada for its financial support during the period of my graduate studies.viiiPart IShip Hull Design OptimizationSummaryThe current study originated from a series of seakeeping experiments on twenty-sevendifferent ship hull designs. The experimenters sought a method of predicting perfor-mance for hull designs different from those included in the study. The results of theseexperiments were analyzed by using linear and nonlinear regression (Delampady et al.,1989), and by minimizing three different loss functions (Delampady et al., 1990).Continuing the work on this problem, we wished to consider alternate methods of anal-ysis. Since the Taguchi method provides a well-established, but somewhat controversialprocedure for design optimization, we decided to evaluate it as a possible alternative.In this report, we describe the initial design problem and the 1989 and 1990 analyses ofthe experimental data. We then provide a detailed description and critical evaluation ofthe Taguchi methodology. Our purpose is to assess its merits in general and to determinewhether it would be suitable for use in the ship hull design problem. We conclude thatthe Taguchi method is not logically sound, and is inappropriate for use in the ship hulldesign problem.11 IntroductionWhen ships are designed, marine architects must consider their performance in all kindsof weather. They seek ship hulls which can be operated in a selected range of speedswithout exhibiting any objectionable tendencies, such as pitching, rolling, heaving, andslamming.New designs are often tested by towing models in tanks at various speeds and underdiverse water conditions. A series of such seakeeping experiments was conducted bythe Institute for Marine Dynamics of the National Research Council of Canada; a fullsummary of the resulting data can be found in Murdey and Butler (1985). The experi-menters sought a method of predicting performance for hull designs different from thoseincluded in the study, based on the results of these experiments.The data from twenty-seven different ship hull designs, of which twenty-five possessedproven satisfactory performance at sea, were analyzed by Delampady, Kan, and Yee(1989, hereafter DKY89) and by Delampady, Gu, Ma, Meloche, and Molyneux (1990,hereafter DGM90). Our purpose in the current study is to investigate an alternativemethod of analysis. Since the Taguchi method is a popular, but controversial procedurefor design optimization, we have chosen to review the methodology to determine whetherit would be an appropriate alternative.In the following sections, we describe the 1989 and 1990 analyses. We then explain themechanics of the Taguchi method and discuss its possible application to the currentproblem. Finally, we examine criticisms of the methodology.22 Initial AnalysesBefore describing the initial analyses, we introduce the following standard naval archi-tecture terms in non-technical language to give an idea of the meaning of the variableswhich will be used (for further details, see Comstock (1967)):Fn, Froude number: V/VTL, where V = speed of ship, g = acceleration due togravity, and L = length of ship; En provides a measure of the ship's speed inrelation to its size,Tz , non-dimensional zero crossing wave period: a measure of wave conditionswhich includes wave period and height,L 2/BT, length-displacement ratio: the ratio of length to displacement where B =beam of ship measured at load waterline, and T = draught, the depth of waterrequired for the ship to float,B/T, beam-draught ratio: the ratio of beam to draught,CB, block coefficient: V/LBT, where V = volume of ship's displacement; CB repre-sents the ratio of the volume of displacement to the volume of a rectangular solidof equal length, breadth, and depth,Cw , waterplane coefficient: (area of waterplane)/LB, where the waterplane is thehorizontal plane through the ship at load waterline; Cw is the ratio of the area ofthe waterplane to the circumscribing rectangle,pitch motion: rotational motion of a ship about its lateral axis,heave motion: motion of a ship along its vertical axis,bow acceleration: a measure of vertical bow acceleration,bow relative motion: a measure of the motion of the ship's bow incorporating verticalacceleration and rotational motion about the ship's longitudinal axis,thrust increase: the increase in forward driving force on the ship.3The first analysis of the experimental data was conducted by DKY89. The explanatoryor predictor variables used by the experimenters were:Fri : Froude number andil : non-dimensional zero crossing wave period,together with the design variablesL2/BT : length-displacement ratio,B IT : beam-draught ratio,CB : block coefficient, andCw : waterplane coefficient.Their response variables were denoted by:yi pitch motion,Y2 heave motion,y3 bow acceleration,y4 bow relative motion, andy5^thrust increase.High levels of all these responses are considered undesirable. They cause discomfort forthe crew, impede the ability to maintain headway, and can coincide with dangerouslyhigh forces on the ship's hull, leading to damage and possible disaster. There may beother ramifications as well, depending on the type of ship. On a naval ship, for example,the limits for pitch might be critical if missiles were being tracked or sonar were used; onan aircraft carrier, the amplitude of vertical motion at the end of the flight deck wouldbe important.DKY89 used both linear and nonlinear regression to relate each of the response variablesseparately to the set of explanatory variables. In the case of bow relative motion, thelinear model provided a slightly better fit; for pitch motion and heave motion, the non-linear model proved to be better. For the other response variables, no real improvementwas provided by the nonlinear model over the linear model.The model for thrust increase was by far the worst fitted model among those for thefive responses considered. It was observed that there was a general lack of pattern in4the variation of thrust increase with respect to the explanatory variables. Furthermore,it was felt that the explanatory variables considered may not have been adequate formodel building.The 1989 study concluded with three recommendations:1. Since the responses were measured simultaneously, they could exhibit an inter-esting multivariate structure, which should be exploited in the construction of anoptimal design.2. Only five of the most important response variables were analyzed; the analysisshould be extended to the rest of the variables which appeared in the experiments.3. Since there are many design variables and environmental factors involved, it wasrecommended that the response surface model of Weerahandi and Zidek (1988) beused.A subsequent study (DGM90) was directed towards finding the ship model with the bestpossible performance. In this study, d was used to denote the vector of hull form designvariables, which again consisted of L2 /BT, B/T, CB, and Cw; a denoted the adjustablevariable Frt and c the uncontrollable variable Pz . Frz was called an "adjustable" variablesince under actual operating conditions the value of Fn, in other words the ship's speed,would be controlled by the captain; Pz was considered an "uncontrollable" variable sincethe operator would have no control over sea conditions.The response variable used by DGM90 was the vector R = (R1, R2, R3) ' where thecomponents of R were:R1 :^heave motion;R2 :^bow acceleration; andR3:^thrust.The set of response variables used here did not include pitch and bow relative motion,since the experimenter believed the information provided by these responses was in-cluded in R. Under the assumption that d, a, and c provide the best explanation ofvariation in R, we have R = R(d, a, c).5The investigators considered a loss function L, where L = L(R(d, a, c)). They sought adesign d which minimized L, bearing in mind that the same design might not be optimalfor all a. Furthermore, since c is a combination of a large number of uncontrollablefactors, they decided to treat c as a random variable, and to average L with respectto the probability distribution of c. The average loss due to the use of a ship hullof design d at speed a was defined to be q(d, a) = Ec[L(R(d, a, C))1, where C is arandom variable which assumes values c, and EC denotes the expectation with respectto the probability distribution of C. The problem of locating the optimal design d nowconsisted of finding d which produced the minimum value of 0. Of course, it might befound that a different optimal design d existed for each speed a at which the ship couldbe operated. This question will be considered shortly.Multiple response, multiple linear regression techniques were used to model R as afunction of c, where the linear model of DKY89 was used. Three different loss functionswere considered. In the notation used by DGM90, these three functions were:3L (Xi ) X21 X3) =^ixil,i=13L(x i , x 2 , x3 )^Exi2, andi=13L (Xi, X2, X3) = Eexpax i p.i=1Performance was measured by averaging 0 over the various Froude number values (de-noted by a), since they were all considered to be of equal importance.The investigators found that the relative ranks of hull models, determined by minimizing0 over d, were not sensitive to the choice of loss function from the three functions givenabove. In addition, they concluded that hull length was not an important factor in thedetermination of the optimal design, because the relative rankings remained essentiallythe same for different hull lengths.The hull forms with good performance all had similar design dimensions. In particular,L2/BT of 150 and BIT of 5.20 seemed to provide the best performance. However, inview of the fact that there would certainly be interaction between these two parameters,it could only be said that the best design possessed parameter values close to these.6The 1990 study recommended that tests be done on hull forms close in dimensionsto those which were found to be optimal. Furthermore, it suggested that in futureanalyses, the hull length be taken as 120 m., and the expected losses determined onlyfor the squared error loss function. A nonlinear model could be used in place of thelinear model used here.Because the focus of this project was the quality of ship hull performance, we felt thatit would be interesting to consider the use of the Taguchi method for the solution of thisdesign problem. The Taguchi method is a well-known methodology for experimentaldesign and analysis; its purpose is to incorporate quality goals in the initial design stageof a product rather than to make costly adjustments to achieve quality at a later stagein the manufacturing process. This emphasis suggested that the method might be veryappropriate for the current situation.73 The Taguchi Method of Parameter DesignIn this section, we describe the basic elements of the Taguchi method of design opti-mization, beginning with the ideas motivating the method. Our description is basedprimarily on readings of Byrne and Taguchi (1987), Ross (1988), and Tunner (1990).3.1 Quality and the Loss FunctionTraditionally, manufacturers and consumers have had rather distinct views of quality.The consumer will be concerned with utility, durability, value, aesthetics, and variousother attributes. Although initially the manufacturer must surely have the same con-cerns, by the time a product reaches the manufacturing stage, the manufacturer willhave translated these criteria into the simple issue of whether the product meets itsdesign specifications.Generally speaking, it is intuitive that quality is improved as variation is reduced. In anengineering environment, this intuitive view can be made precise. Parts designed to fittogether must each be the correct size or the final product will work badly or perhapsnot at all. Engineers are also aware, however, that it is impossible to reduce variabilityto zero. For this reason, engineers specify both a target value and acceptable limitswhen establishing specifications.If the target for a dimension were 1000 mm., for example, the design specifications mightthen be 1000 ± 2 mm.. In other words, the lower specification limit (LSL) is 998 mm.and the upper specification limit (USL) is 1002 mm.. This implies that a part which is997.9 mm. long (part A) is completely unacceptable, while a part which is 998.1 mm.long (part B) is exactly as good as one with the precise target length of 1000 mm(part C). Furthermore, it is implied that the production of part B will cause no loss tothe company as a result of its deviation from the target value.From the customer's point of view, there is very little difference between part A andpart B, since the customer is concerned with the final function of the part. There might,however, be a good deal of functional difference between part B and part C.8Loss LossNo LossA B^CLSL USLTargetFigure 1: Traditional loss function (Ross, p. 42)To the manufacturer using the artificial LSL, on the other hand, there seems to be anenormous difference between parts A and B and no difference at all between parts Band C. The manufacturer will incur scrap or rework costs for part A if inspection is onehundred per cent effective, or warranty costs if inspection is not completely effective.Part B will supposedly generate no costs whatsoever, and part C will definitely generateno costs.The "traditional" loss function reflects this viewpoint. It is shaped like a "goal-post": itsvalue is zero inside the product specifications and takes a "step" jump at the specificationlimits, as can be seen from Figure 1.In practice, of course, this loss function represents a gross oversimplification. The sizeof a car door, for example, may be within specifications but not quite on target. Theloss to the consumer could range from simple annoyance at having a door which fitssomewhat imperfectly to the inconvenience of water leaks. The resulting loss to themanufacturer could be damage to its reputation and possible loss of market share.It is therefore necessary to reconcile the customer's and manufacturer's views by usinga loss function which reflects both perspectives. It is also important to note that ascustomers become more demanding, they increase pressure to reduce variability.A simple parabolic curve, as shown in Figure 2, has an intuitive appeal, especially whenspecifications are symmetric about a target value. The loss is very small close to the9j Increasing LossFigure 2: Parabolic loss function (Ross, p. 42)target and grows at an increasing rate as one moves away from the target. The curveflattens to a value equivalent to scrap cost outside the specifications since it is reasonableto assume that generally speaking, products which fall outside the specification limitswill be scrapped.3.2 The Taguchi ApproachGenichi Taguchi is a Japanese engineer whose ideas have had a significant effect onJapanese and more recently North American technology, beginning with his work in theearly 1950's. To incorporate his idea that loss increases with increasing distance fromthe target value, Taguchi uses a quadratic loss function like that shown in Figure 2:L = K (S2 + (X — T) 2 ),whereL =K =S2 =X =T =average loss per unit produced,cost constant,process variance,average specification value, andtarget specification value.10Furthermore, if one definesc = loss associated with a unit of product produced at a specification limit(assuming that the loss for a unit at target is zero), andd = distance from the target to the specification limit,then,K = c/d2 (Tunner, 1990, p. 58).While it may be difficult to produce a valid estimate of c, it should be possible for expertsin the subject area to derive a reasonable value. One group of engineers has suggestedthe rule of thumb that c be one-tenth of the selling price of the item in question. (seeTunner, 1990, p. 59)The quadratic loss function underlies Taguchi's entire approach to quality engineering.The function shows the economic advantage of reducing variation while incorporatingcustomer requirements. Reducing variation and staying close to the target are bothimportant, since either or both of these achievements will lead to reduction in total loss.In fact, it is often the case that cost will be reduced at the same time that quality isimproved.Furthermore, Taguchi divides system quality engineering into two parts: off-line and on-line quality control. He stresses the need to engineer quality into a product in the designstage in order to minimize the costs of maintaining quality during the manufacturingprocess. The off-line process will establish a design target for minimum loss and reducevariance by design before production actually starts. Later, the on-line process will holdthe average close to the target.According to this scenario there are three steps in the (off-line) engineering optimizationof a process: system design, parameter design, and tolerance design. System designinvolves the selection of materials, parts, equipment and tentative parameter values.The parameter design stage determines the optimal levels of parameter values; values forproduct parameters and process factors are chosen which are least sensitive to changesin noise factors. This is the crucial stage for attaining high quality at minimal cost.If parameter design does not reduce variation sufficiently, tolerance design must be11used. This usually implies greater cost, since it means buying better materials, parts,or machinery in order to tighten tolerances on product parameters or process factors.To clarify these ideas, consider the following simplified example. In the system designstage for the manufacture of ceramic tiles, the type of kiln and raw materials would beselected, along with tentative values for proportions of raw materials. Suppose it wereknown that the temperature within the kiln varied by a specific amount; this variationwould be considered a noise factor. In the parameter design phase, proportions of rawmaterials would therefore be selected which would result in minimal sensitivity to thistemperature variation while producing tiles which met the desired specifications. If theparameter design phase were omitted, it might later be necessary to make expensive kilnmodifications (tolerance design) in order to reduce kiln temperature variation, therebyachieving the necessary tile quality.Note that even in cases where parameter design has been omitted initially, it is possibleto go back to this stage after production has begun. In this way, one can still hope toavoid incurring the expense of tolerance design. The tile production example used abovewas taken from a real case in Japan in 1953 (Byrne Si Taguchi, 1987, hereafter BT87,p. 20). By choosing better factor levels, the percentage of size defects in the ceramictiles produced was reduced from 30% to less than 1%.In North America, engineers commonly jump from system design to tolerance designsince they are accustomed to spending money to obtain the required performance levels.In addition, Taguchi methods are often applied to a current production problem toimprove quality, when they could have been used much earlier in the design phase.Taguchi's idea of reducing variation is further refined by the emphasis on doing thisin the most cost effective manner possible. The important factors in a process areseparated into control factors and noise factors. Control factors are those factors which amanufacturer can control relatively easily. Noise factors are those which a manufacturercannot or does not wish to control, perhaps due to the high cost of such control. Thegoal of the parameter design process is to reduce variation by choosing levels of thecontrol factors which will make the response less sensitive to changes in noise factorlevels, since noise factors, by their very nature, are usually more expensive to control.Control factor levels are chosen by conducting appropriate experiments.125 7 5 10Change in HeightConsistency Batch-to-Batch3020Change inHeight(growth)(0.0001 in.)10AmmoniaMeasuringPrecisionAmount of Ammonia (ft 3 )Figure 3: Result of Heat Treatment Experiment (Ross, p. 47)Figure 3 illustrates the result of a designed experiment which was used to improvea heat treatment process. The goal was to meet the height requirements for a camfollower. The precision of the ammonia measuring system, i.e., variation in the amountof ammonia, was treated as a noise factor, since the cost of improving it would havebeen significant. The experiment indicated that the amount of ammonia had a stronginfluence on controlling the change in height experienced in the process. The previousproduction standard of 5 cubic feet of ammonia resulted in a large variation in heightfor the cam followers. When the amount of ammonia was increased to 8.75 cubic feet,the variation in height was reduced substantially. In other words, the process wasless sensitive to variation in quantities of ammonia when the amount of ammonia wasincreased. The cost of the additional ammonia was very small. Here again, the conceptof using parameter design to find parameter levels (in this case amount of ammonia)which were relatively insensitive to noise (variation in amount of ammonia) resulted inan improved process at relatively low cost.Note that "the search for interactions among controllable factors is de-emphasized, al-though there are exceptions. The key to achieving robustness against noise is to discoverthe interactions between controllable factors and noise factors. Specific interactions be-tween controllable factors and noise factors need not even be identified. As long as13Controllable Factors LevelsA. Interference Low Medium HighB. Connector Wall Thickness Thin Medium ThickC. Insertion Depth Shallow Medium DeepD. Percent Adhesive inConnector Pre-dip Low Medium HighUncontrollable Factors LevelsE. Conditioning Time 24h 120hF. Conditioning Temperature 72°F 1500FG.^Conditioning Relative Humidity 25% 750/oFigure 4: Factor levels for elastomeric connector experiment (BT87, p. 21)the noise factors are changed in a balanced fashion during experimentation, preferredparameter values can be determined using an appropriate signal-to-noise ratio." (BT87,p. 21) The treatment of interactions and the use of signal-to-noise ratios for analysiswill be discussed later. First, we will consider the design of experiments.3.3 Experimental DesignThe Taguchi method uses a crossed array design which is most easily explained by theuse of an example. In this example (BT87, pp. 21-22), researchers aimed to attach anelastomeric connector to a nylon tube so as to maximize pull-off force, the force requiredto pull the connector off the tube. The researchers selected four important controllablefactors: A) interference, B) wall thickness, C) insertion depth, and D) percent adhesive,and tested each at three levels. They also selected three noise factors: E) conditioningtime, F) conditioning temperature, and G) conditioning relative humidity, each to betested at two levels. Figure 4 shows the factors and levels. As noted earlier, noise factorsare generally uncontrollable or too expensive to control during normal operations, butcan often, as here, be controlled for experimental purposes.The crossed array is formed from an inner array and an outer array. The inner array isan orthogonal array containing the levels of the controllable factors to be tested. An L9array (shown in Table 1) was selected for this purpose as the most efficient orthogonal14FactorNo.ABCD1 1 1 1 12 1 2 2 23 1 3 3 34 2 1 2 35 2 2 3 16 2 3 1 27 3 1 3 28 3 2 1 39 3 3 2 1Table 1: L9 array for elastomeric connector experiment (BT87, p. 22)design for four factors at three levels. Its task was to find the best of the 3 4 = 81 possiblefactor combinations.The outer array is an L8 orthogonal array containing the levels of the noise factors tobe tested (see Table 2). The columns labelled E, F, and G show the levels at whichthe noise factors are to be set. The other columns are used in the analysis to estimateinteraction effects and experimental error. "The most important use of this array is todeliberately create noise to identify the controllable factor levels that are least sensitiveto it. Usually the effects of interactions among the noise factors are not estimated, butthe experimenters in this case allowed for such interactions, as they felt the informationmay be valuable." (BT87, p. 22)The two arrays are then combined into a crossed array (see Figure 5). Each controlfactor combination is run with each noise factor combination, resulting in a total of 72runs.In this case, it was considered reasonable to run 72 experiments. Had it been desirableto reduce the number of runs to 36, an L4 outer array could have been used insteadof the L8 array, at the expense of being able to examine interactions among noisefactors. The number of runs could have been reduced even further to a total of 18, by158 7 6 5 4 3 22 2 2 22 2 2 21 2 2 2 22 1 2 2 1 2 11 2 2 2 1 22 2 1 2 22 1 2 2 2120h 120h 120h 120h 24h 24h 24h 24h150E 150F 72F 72F 150F 150F 72F 72F75% 259 75% 25% 75% 25% 75% 25%1^.128.6^22.6^22.7^23.1^17.3^19.3^19.9^16.1J< C/)ZZ 0Ezx 9 3 3 2 1^High^Thick^Medium^LowL9 ARRAY  I Inver - Wall Ins. Percentference Thickness', Depth AdhesiveHigh Me on S a ow High(A^ (C)^D1Cond.Temp.(F1Cond.S.M.Cond.Time(F.)ExC(G)ExFFxGeEL8 ARRAYS/NRatio(db)26.152No.FactorE FExF GExGFxG e1 1 1 1 1 1 1 12 1 1 1 2 2 2 23 1 2 2 1 1 2 24 1 2 2 2 2 1 15 2 1 2 1 2 1 26 2 1 2 2 1 2 17 2 2 1 1 2 2 18 2 2 1 2 1 1 2Table 2: L8 array for elastomeric connector experiment (BT87, p. 22)Figure 5: Crossed array for elastomeric connector experiment (BT87, p. 22)16using an outer array which simply compounded the noise factors to two noise levels, a"nominal" and a "worst" condition. In addition to the restriction produced by the 36run design, this design would also make it impossible to examine specific interactionsbetween controllable factors and noise factors. Their identification is not generallyconsidered necessary, however, due to the method of analysis used. In the followingsection we continue our exploration of the Taguchi method with a description of theanalysis procedure.3.4 AnalysisTaguchi's standard method of analysis examines a signal-to-noise (S/N) ratio, wherethe signal is considered to be the mean and the noise is the standard deviation. Thereare many different S/N formulae, but there are three which are generally used, one foreach of three different response classifications: "the higher the better", "the lower thebetter", and "nominal is best". In each of these situations, the S/N ratio is constructedin such a way that the S/N ratio is to be maximized. The idea is that each of theseratios offers the best balance between optimizing the mean response and minimizing theresponse variation for the circumstances in which it is intended to be used.Continuing with the elastomeric connector example, pull-off force falls into the "higherthe better" classification, so the S/N ratio used is:S/N (in dB) = —10 log [1 (1^17., —yi2 + --y22 + • • • + —y1.2)] •where Yi, Y2, • • • , yn refer to the n observations (in this case n = 8) for each setting ofthe controllable factors.The analysis can be done in several ways. One common technique is to use F tests ina statistical analysis of variance (ANOVA) to find the statistically significant factors.Taguchi, however, often teaches engineers to use a conceptual approach which involvesgraphing instead. For each controllable factor, the average S/N ratio for each level isplotted on a graph (see Figure 6). The graphs show the factors which seem to have thestrongest effect and the optimum setting for each factor when the mean and varianceare considered together. Mean response plots (as shown in Figure 7) can be examined as1727.0 -26.0 -25.0-24.0 -IThin Med Thick(B) Wall Thickness27.0-26.0-25.0-24.0 -27.0 -26.0 -dB^//^`25.0 - f24.0 -Low^Med^High(A) Interference27.0 -26.0 -25.0 -24.0 -Shallow Mel d Deep^Low Med High(C) Insertion Depth^(D) Percent AdhesiveFigure 6: Average S/N ratio for each control factor level (BT87, p. 23)well. These can supply additional information, since they can show whether the meanor the variance has the stronger effect on the S/N ratio.If a factor does not have a strong effect, it is less important to choose the optimum settingfor that factor from the graphical analysis, since there is apparently little differencein the response produced at each level. Instead, other considerations such as cost orconvenience may be used in the selection process.Although it is claimed that performing the S/N analysis eliminates the need to examinethe specific interactions between controllable factors and noise factors, this is sometimesdone graphically as well, in order to obtain more information about the problem. Inthe example in question, this graphical analysis was done for all of the control x noiseinteractions. The results did not alter the decisions made through the S/N and meanresponse analyses.A separate study had been conducted on ease of assembly. When the results of the twostudies were combined (which involved some compromises), the optimal control factor1821.0 -20.0 -19.0 -18.0 -21.0 -20.0-19.0-18.0 -Low Med^ Id Hl jgh^Low^Med High(A) Interference^(B) Wall Thickness21.0 -20.0 -19.0 -18.0-Shallow Med Deep21.0-20.0-19.0-18.0 -LOw Med High(C) Insertion Depth^(D) Percent AdhesiveFigure 7: Mean responses for each control factor level (BT87, p. 23)settings were determined to be A2, B1, C2, and D1 . This control factor combinationhad never been run during the experiment, so a five-run confirmation experiment wasconducted, which did indeed produce increased pull-off force.If researchers feel that for some reason none of the three S/N ratios given above isappropriate for their situation, a common alternative divides the analysis into two steps:1. Examine the S/N ratio, S/N (in dB) = 20 log(y/s), where y and s are the meanand standard deviation of the response from a given experimental condition. Thefactor levels which correspond to the highest S/N ratio are chosen to minimizevariation.2. Factors which do not appear to be significant in the S/N analysis do not havea strong influence on variation, and can therefore be chosen to adjust the mean,either towards a target level, or towards a maximum or minimum.This approach puts more weight on variation, whereas the standard S/N ratios put moreemphasis on the mean.194 The Taguchi Method and the Ship Hull DesignProblemWe shall now describe the author's interpretation of how the Taguchi method mightbe applied to the ship hull design problem; this is done to ensure that the method is,in fact, applicable to this problem. After our critical review of the Taguchi method inSection 5, we shall determine whether we consider its use in the ship hull design problemdesirable.Recall that in the 1990 study of the ship hull problem, the design variables were:1,2 I BT : length-displacement ratio,B IT : beam-draught ratio,CB : block coefficient, andCw : waterplane coefficient.In Taguchi parlance, we would call these the control factors. The set of explanatory orpredictor variables comprised these control factors together with:Fa : Froude number, and11, : non-dimensional zero crossing wave period.We would call t a compound noise factor, but there might be some uncertainty aboutthe classification of Fn, the Froude number. Since this factor is adjusted when the shipis in operation, it is certainly not a noise factor. On the other hand, it should not beclassified as a control factor because that would result in the choice of only one levelof Froude number when the final analysis is complete. The fact is that we seek goodperformance at all Froude number values. A workable solution might be to conductseparate full experiments for each of various Froude numbers and then average theresults, assuming that good performance at each of these speeds is equally important.This is a similar approach to that taken in the second study discussed earlier (DGM90).In that study, it was found that little information was lost in this averaging process.Thus far we have dealt with the control and noise factors. Next we turn to the re-sponse variable, which must be scalar valued for the application of the Taguchi method.20However, the response in the DGM90 study was vector valued R. (R1, R2, R3)' withcomponents:R1 : heave motion;R2 : bow acceleration; andR3: thrust.We must therefore modify our response in some way to make Taguchi's method appli-cable. One could optimize each of the three responses separately, of course, but doingthat would have the same shortcomings as it did in the previous analysis.We would prefer to construct a function for an overall performance measure which wouldinclude each of the components of R, and possibly other performance measures as well.This would need to be done by a specialist with extensive knowledge of marine archi-tecture. Historically, marine architects have not developed such a function; productionof new designs has relied more on experience than on quantitative tools. Recently, how-ever, there has been some effort to develop a generally accepted performance index (seeLloyd and Andrew, 1977). It would be interesting to conduct a study with such anindex as the response to be optimized. Without such an index, we have not been ableto proceed further at this time. However, this might well be the subject of future work.Once the response is defined, one could construct the inner and outer arrays accordingto a Taguchi orthogonal design. The outer array would be very simple, since thereis only one compound noise factor. The analysis could be done using the S/N ratiofor the "highest the best" scenario. There is a potential problem which may arise ifthe optimization analysis were conducted with just the currently available data. Inthe Taguchi method, one envisions the construction of orthogonal arrays and a seriesof experiments on those arrays. But for the available data, the settings at which theexperiments were run may not fit into the envisioned structure. This would need to bechecked.We have learned the mechanics of the Taguchi method and discussed fitting our probleminto the Taguchi framework. In the next section, we will examine the method criticallyto determine its potential value in the solution of our problem.215 Criticisms of the Taguchi MethodThe motivating philosophy behind the Taguchi methodology appears sound; in qualitycontrol literature it is widely applauded (see Ryan, 1988 or Nair, 1992, hereafter N92,for example). The same cannot be said of the methodology developed to implement thephilosophy, however. This area is rife with controversy.The main criticisms of the method, as found in the author's review of subject arealiterature, seem to be:1. Experimental Design:(a) A crossed orthogonal array design is seldom the most efficient experimentaldesign. It often needs more experimental runs (and thus incurs higher costs)than are really necessary.(b) Ignoring control factor interactions as a matter of course is not a good policy.(c) The design does not allow for sequential experimentation.2. Analysis(a) Marginal Average Plots:(i) The examination of plots of marginal averages is not a sound techniquefor determining the best control factor settings.(b) Response Function:(i) The signal-to-noise ratio is not the best response to examine, and canjust complicate the situation unnecessarily.(ii) The loss function should be used directly in the optimization processrather than disappearing from the situation after motivating the method.3. General Approach:(a) The methodology suffers from a preoccupation with optimization rather thanseeking to accomplish optimization through an understanding of the physicalprocess.We now examine each of these criticisms more closely.225.1 Experimental Design5.1.1 EconomyFor each combination of control factors in the inner array, observations are made at allcombinations of the noise factor settings in the outer array. If In and n are the respectivesizes of the inner and outer arrays, the product array will have mn runs. If control andnoise factors are included in a single array, however, the total number of runs requiredcan be significantly smaller. Several authors, including Shoemaker, Tsui, and Wu (seeN92), have suggested the use of a combined array for reasons of economy.Furthermore, we must recognize that the distinction between control and noise factorsis somewhat arbitrary. As far as the physical system itself is concerned, this distinctionis merely an artifice. It makes more sense to use a design which simply includes allthe appropriate factors in the most efficient manner possible. According to Lorenzen(N92, p. 150), the combined array usually has a better confounding pattern than theproduct array, as well as superior robustness properties. The combined array, however,may sometimes require more control factor combinations; if control factor settings aremore expensive to change than noise factor settings, the combined array can be moreexpensive to use in these cases.5.1.2 Interaction EffectsSacks and Welch (N92, p. 140) suggest that Taguchi has ignored interactions to minimizethe number of runs needed by a product array. This seems indefensible and has beenwidely criticized by other authors as well, such as Ryan (1988), Hendrix (1991), andBox (N92), to name a few.Byrne & Taguchi (BT87, p. 24) claim that "performing the S/N analysis generallyeliminates the need for investigating the specific interactions between controllable factorsand noise factors," although they offer no argument to support this statement. It isdifficult to see the logical connection. In fact, they go on to say that the investigationof interactions may provide additional insights into a problem. As George Box has23noted (N92, p. 146), "because engineers have traditionally relied on one-factor-at-a-time experimentation, main effects will often have already been put to use, and it willbe the unexpected interaction that is waiting to be discovered and sometimes to beexploited with dramatic results." It seems most peculiar to maintain that we mustexperiment to find the factors with significant main effects, while claiming at the sametime that we know the control factor interactions are unimportant. Furthermore, wesearch for interactions between control and noise factors, while ignoring interactionsbetween control factors. Since control and noise factors are alike in a mathematicalsense, the absence of control factor interactions seems implausible; certainly they cannotbe expunged by merely classifying certain factors as control factors.Even if we accept the use of product arrays, we can often improve Taguchi's inner arraydesigns without incurring the cost of additional runs. Taguchi advocates the use oforthogonal arrays for the control factors and has published a number of these arraydesigns, some of which are suboptimal. For example, Ryan (1988, p. 35) describes onecase where Taguchi and Wu show an L8 array for investigating four main effects and twotwo-factor interactions, where three of the main effects are confounded with two-factorinteractions. Instead of using this array one could simply use a 24- i fractional factorialwhere main effects are confounded with three-factor interactions. This is a far betterdesign, with exactly the same number of experimental runs.Strictly speaking, here Ryan's criticism is of a specific implementation of the Taguchimethod, rather than the method in general. It does, however, illustrate that it is veryimportant to be aware of the aliasing structure of an array. If a main effect is confoundedwith a two-factor interaction and the experimenter is unaware of this, it would be easyto conclude that the main effect is significant when in fact the interaction effect issignificant, and vice versa. Section 5.1.3 will explore this idea further.5.1.3 Sequential ExperimentationThe Taguchi methodology requires that the experimental design be chosen, run, and thenanalyzed to determine the optimal control factor settings. When a fractional design isused in a "one-shot" approach, however, one must be very careful. It is important to24make a correct initial decision about the effects to be tested. If this decision is incorrect,everything which follows may be completely misleading due to aliasing. Unless we use afull factorial design, it seems sensible to take an iterative approach to "zero in" on theimportant factors. It is therefore often useful to do a preliminary screening (fractionalfactorial) experiment based on the researchers' initial views. The results can dictatewhich fraction should be run next in order to continue the investigation without havingimportant effects aliased with each other.Barad, Bezalel, and Goldstein (1989) give an excellent example of a sequential exper-iment. It stresses the implications of aliasing as well as the usefulness of preliminaryresearch to determine which variables are important. In the example given, productdevelopment engineers tested six battery design parameters at two levels each in orderto find out which of them significantly affected battery capacity. A full factorial experi-ment would have required 64 runs, but since the engineers felt that the most importanteffects were the six main effects as well as the BD, CD, and EF interactions, they wereable to do a quarter-fraction screening experiment consisting of 16 runs.The analysis of the quarter-fraction indicated that the BD interaction was very signifi-cant; this was surprising, since neither of the main effects B and D was significant. Thealias structure showed that BD was aliased with AE, which had not been consideredimportant in the initial discussion. Further investigation found that the AE interac-tion and several others which had been considered unimportant earlier were, in fact,significant.5.2 Analysis5.2.1 Marginal Average PlotsIn general, plotting marginal averages to determine control factor settings simply willnot work. This practice is comparable to conducting an experiment where the factors arevaried one at a time. Plotting marginal averages is not useful unless all the interactioneffects are zero, and as seen earlier, it is completely unreasonable to assume that thiswill necessarily be the case. Ryan (1988, pp. 34-35) gives a simple example for an2518141061410168B HIGHB LOWINTERACTION PROFILE^ MARGINAL AVERAGESResponse^ Average^ AverageVariable Response ResponseLow^High A Low High A Low^High BFigure 8: Interaction and marginal average plots for the 2 2 design (Ryan, pp. 34-35)Average^ AverageResponse Response1110138Low^High ALow^HighBFigure 9: Marginal average plots for the revised 2 2 data (Ryan, p. 35)unreplicated 22 design. Figure 8 shows the interaction and marginal average plots forthe data. The lines on the interaction plot are parallel, indicating that there is nointeraction between the two factors. The plots of marginal averages correctly indicatethe choice for the highest response: AhighBhigh = 18.Now suppose the response at AhighBhigh had been 12 instead of 18. The new marginalaverage plots are shown in Figure 9 and would again suggest the choice of AhighBhigh .In this case, there is an interaction between the two factors and the marginal plots donot give the correct solution, since the maximum response is in fact AhighBlow = 14.The inescapable conclusion is that if we cannot rely on marginal average plots to giveus the correct answer when all factor level combinations are included in an experiment,26we certainly cannot expect them to do so in a fractional design situation.5.2.2 Response FunctionSignal-to-Noise Ratio To illustrate the inadequacy of the signal-to-noise ratio, Hen-drix (1991, p. 75) re-examines the example cited earlier by Phillip Ross (1988, see Fig-ure 3). In that example, high levels of growth were desired. It happened that high levelsof ammonia maximized growth (the signal, in Taguchi's terms); in addition, high levelsof ammonia minimized growth variation (the equivalent to Taguchi's noise). Since highammonia levels maximized signal and minimized noise, they would, of course, maximizethe signal-to-noise ratio. Note that in this example, response optimization coincideswith low sensitivity to variation in the controllable factor.Suppose, on the other hand, that low growth and low noise were favoured. Since lowgrowth coincides with low ammonia levels and low noise with high ammonia levels, somesort of compromise is necessary. The balance depends the importance placed on eachof these two objectives, but the solution will certainly be some intermediate ammonialevel.Taguchi gives several different signal-to-noise ratios for use when the response is tobe minimized, but these alternatives are each simply arbitrary ways of weighting thetwo goals. It would be much simpler to look at the actual graph itself and determineone's own weighting scheme rather than to look at S/N ratios (which in this case areresponse/slope ratios). The situation becomes even more complicated when there areseveral controllable variables.Figure 10 shows the response plots for an example with two controllable variables, X1and X2, where the goal is to minimize the signal Y. Low levels of X 1 are preferable sincethey minimize Y, but the decision for X2 requires some sort of compromise since highlevels of X2 minimize noise but low levels of X2 minimize the signal. Re-expressing thedata as an S/N ratio will not change the fact that a perfect solution cannot be found.Certainly it will help to cloud the picture. Here again, it is simpler just to examine thegraph.27X2 RigllX1Figure 10: A simple example with two controllable factors (Hendrix, p. 76)Hendrix minces no words in his assessment of the utility of the signal-to-noise ratio(p. 76):For those situations in which the response (signal) and the slope (noise) movein agreement with the objective ... , the Taguchi approach embodied as asignal-to-noise ratio will work — but so would a simple graphical presenta-tion of the data. Reducing the data to a signal-to-noise ratio is unnecessary,because it doesn't contribute anything and only confuses the issue.For those situations in which the response (signal) and the slope (noise)are antagonistic ... , reformatting the data can only lead to an arbitrarycompromise. The nature of this compromise is obscured by formatting thedata as a signal-to-noise ratio. Disguising the slope (rate of change) as astandard deviation merely adds to the confusion.At best, one might say that the various S/N ratios used appear to select an arbitrarybalance between the mean and variance which is appropriate for a reasonable numberof cases; in the process, what is really going on is hidden from view and the analyst'soptions are restricted unnecessarily. Not only does the S/N ratio "throw away" someof the information provided by the data, but it complicates the situation needlessly. Ifwe model the quality characteristic of interest rather than the S/N ratio, we will havemore background knowledge about the system to start with, and will gain more insight28into how the factors affect the response. In addition, the quality characteristic is ofteneasier to model than an S/N ratio, resulting in a better model which leads to a morereliable optimization and ultimately, a better design (Sacks and Welch in the article ofN92, p. 147).Furthermore, although compromises between response level (signal) and sensitivity(noise) are addressed by the Taguchi approach, the problem of having to achieve a com-promise among several responses, although very common in practice, is not addressedin Taguchi literature.Loss Function Taguchi's ideas about the loss function and its proper form seemsound. Once the loss function has been postulated, however, it seems odd not to useit directly in the experimental analysis. Instead, its role is transient; it is used onlyto motivate the general approach. Madhav Phadke (N92, pp. 140-141) relates the S/Nratio to the loss function for the "higher the better" case, and claims that a direct min-imization of the loss function itself would not be very effective in minimizing sensitivityto noise factors. He also states that when the loss function is minimized directly, "thereis an increased risk of interaction among the control factors ... ". One would think thatthe existence of interactions would be determined by the physical situation and not bythe choice of function to be optimized.Even if the S/N ratio is related to the quadratic loss function, however, it must be notedthat the simple quadratic loss function is not always the correct function to use (Tribusand Szonyi, 1989, p. 50). For example, when deciding on the proper length of rope tospan a river, it is much more expensive to choose a rope which is too short than onewhich is too long.As observed previously, another shortcoming of the use of the S/N ratio is that it doesnot allow for multiple responses. These could be incorporated into a loss function whichwas to be optimized directly; for example, the loss functions for the individual responses,appropriately weighted, could simply be added together.To sum up, one could say that it makes more sense to work with the item of interest itself(i.e., either the loss function or the quality characteristic) than an artificial construct.29It appears that the S/N ratio was introduced to simplify the solution. The problem isthat in doing so, it oversimplifies the situation, leading potentially to erroneous results.It may also encourage practitioners to lose sight of their original purpose.5.3 General ApproachThe Taguchi method very nicely lays out a path to be followed in the quest for opti-mization. Unfortunately, it is difficult to chart a course at the outset, when one doesnot know the true destination. As George Box (N92) has pointed out, it is best totry to increase understanding of the process, rather than to try to optimize it withoutinvestigating the underlying structure. One is then more likely to find a true optimumrather than what might simply be a local optimum, and the cost of doing this is notnecessarily any greater. As Box has eloquently stated (N92, p. 131):It seems that Taguchi's experimental strategy is intended only to pick the"optimum" factor combination from a one-shot experiment. Although theimmediate objective may be this, the ultimate goal must surely be to betterunderstand the engineering system. For example, appropriate designs canprovide estimates of those specific interactions between environmental anddesign factors that cause lack of robustness. Once the engineer knows whichthese are and what they can do, he can employ his engineering know-howto suggest ways of compensating for them, eliminating them, or reducingthem. Thus I profoundly disagree with the implications of Shin Taguchi'sclaim that the engineer does not need to "discover the causal relationshipsand to understand the mechanics of how things happen." To believe this isto discount the way the engineer thinks and the whole basis of his educationand experience.306 ConclusionsNot surprisingly, we must conclude that although the Taguchi method may be usefulin some cases, it is, in general, not a sound optimization method. As delineated inSection 5, we have serious misgivings about both the experimental approach and theanalytical procedure.Even if we were convinced that the analytical steps were reasonable and attempted toapply them to the data which were collected in the ship hull experiment, we would firstneed to ensure that the data would fit into a Taguchi-style crossed array. Next, since weare interested in simultaneously optimizing a number of responses, we would be requiredto construct some sort of overall response function which could be forced into one of thethree standard analysis situations. This would certainly lead to many of the attendantproblems described in Section 5.2. Although it would certainly be interesting to seethe results of the Taguchi method applied to the ship hull design problem, we favourmore conventional alternatives, such as response surface methodology, for the purposeof finding a "true" optimum.31References for Part IBarad, M., Bezalel, C., and Goldstein, J.R.(1989). Prior Research is the Key to Frac-tional Factorial Design. Quality Progress, 22, 3, 71-75.Byrne, D.M., and Taguchi, S. (1987). The Taguchi Approach to Parameter Design.Quality Progress, 20, 12, 19-26.Delampady, M., Kan, L., and Yee, I. (1989). Analysis of Seakeeping Data. Technicalreport submitted to the National Research Council of Canada.Delampady, M., Gu, C., Ma, P., Meloche, J., and Molyneux, D. (1990). Optimal HullForms From Analysis of Seakeeping Data. Technical report submitted to theNational Research Council of Canada.Hendrix, C.D. (1991). Signal-to-Noise Ratios: A Wasted Effort. Quality Progress, 24,7, 75-76.Lloyd, A.R.J.M., and Andrew, R.N. (1977). Criteria for Ship Speed in Rough Weather.In Proceedings of the 18th general meeting of American Towing Tank Confer-ence, Volume 2, Cavitation Session, Seakeeping Session, eds. B. Johnson and B.Nehrling, 541-565.Murdey, D., and Butler, J. (1985). Hull form series for fast surface ships: summary ofseakeeping data. Institute for Marine Dynamics, LTR-SH-393.Nair, V.N., ed. (1992). Taguchi's Parameter Design: A Panel Discussion. Technomet-rics, 34, 2, 127-161.Ross, P.J. (1988). The Role of Taguchi Methods and Design of Experiments in QFD.Quality Progress, 21, 6, 41-47.Ryan, T.P. (1988). Taguchi's Approach to Experimental Design: Some Concerns. Qual-ity Progress, 21, 5, 34-36.Taguchi, G., Elsayed, E., and Hsiang, T. (1989). Quality Engineering in ProductionSystems. New York: McGraw-Hill.Tribus, M., and Szonyi, G. (1989). An Alternative View of the Taguchi Approach.Quality Progress, 22, 5, 46-52.Tunner, J. (1990). Is an Out-of-Spec Product Really Out of Spec? Quality Progress,23, 12, 57-59.32Weerahandi, S., and Zidek, J. (1988). Bayesian nonparametric smoothers. CanadianJournal of Statistics, 16, 61-74.33Part IIEndogamous Group ComparisonsSummaryThis report developed from a case which came to the Statistical Consulting and ResearchLab (SCARL) in the Department of Statistics at The University of British ColumbiaHaving lost all her original data, our client wished to know if it were possible to recon-struct her data or to perform a statistical analysis, using the tables of summary statisticsin her possession.She was interested in modelling several different pulmonary functions using various de-mographic and anthropometric measurements as explanatory variables. For each of thedifferent model variants under consideration, her goal was to see if the same parameterswould be valid for all endogamous groups', or if some groups were sufficiently differentthat the parameter values would be different for these groups.Dr. Malcolm Greig of SCARL and I were able to develop an innovative method forthe solution of this problem; the derivation and application of this methodology aredescribed in this paper.'social or tribal groups which do not permit marriage outside the group347 IntroductionThe current study began with a telephone enquiry from a potential client of the Sta-tistical Consulting and Research Lab (SCARL) in the Department of Statistics at TheUniversity of British Columbia She had lost all her original data through a series ofmisfortunes and wished to know if some sort of statistical analysis or reconstructionof the data were possible, since she had managed to retain certain tables of summarystatistics of the data.The client, Dr. Lakshmi Reddy, was interested in modelling pulmonary functions usingvarious demographic and anthropometric measurements as explanatory variables. Foreach of a number of different model variants, her goal was to see if the same parameterswould be valid for all endogamous groups, or if some groups were sufficiently differentthat the parameter values would be different for these groups.In a seemingly hopeless situation, where most prospective analysts had already givenup, Dr. Malcolm Greig of SCARL suggested a route to an innovative solution whichthe author has worked out in detail below. By making certain "working" assumptionsand using a plausible but restricted model, a methodology for testing overall grouphomogeneity in the absence of the original test data is developed. The procedure usesonly the information contained in tables of means, standard deviations, correlations,and group sizes. The author subsequently applied the method to our client's problemwith rather interesting results, which are reported in the following sections.358 Background8.1 The Estimation of Physical FitnessDr. Reddy wished to determine whether the same parameters would be valid for dif-ferent endogamous groups in models of pulmonary functions as a function of variousdemographic and anthropometric measurements. The motivation of her study was thewell-established use of vital capacity estimates as a measure of physical fitness.It is common in industry for job applicants to be given a physical examination andhired only if they are judged to be fit. Vital capacity is a common measure of fitness,but the equipment to measure vital capacity directly is often unavailable. As a result,it is routine for physicians to predict vital capacity based on physical measurementssuch as height and age, using charts designed for this purpose. The charts currentlyused in countries with predominantly Caucasian populations are, of course, based onCaucasians. Their application could result in the rejection of a fit non-Caucasian asunfit, simply because a fit Caucasian of the same height and age might have a largervital capacity. In fact, it is recognized that for orientals, the Caucasian-based standardsshould be multiplied by 0.80. By using other variables in addition to, or instead of, ageand height, Dr. Reddy wanted to develop a system for predicting vital capacity whichcould be applied universally.As part of the research program for her Ph.D. in biological anthropology, Dr. Reddyhad gathered data on a total of 1301 people in India from six different caste groups(since all caste groups in India are endogamous). These included workers from a mopedfactory near the town of Tirupapi in the state of Andhra Pradesh, as well as peoplefrom villages in the immediate countryside. The subjects represented a cross-sectionof people from various backgrounds, and included labourers, housewives, office workers,and professionals. The data were collected over a period of six months in 1980-81 byDr. Reddy and two co-workers, a pulmonary physiologist and an anthropologist. Dr.Reddy intended to use these data to investigate her idea after the completion of herPh.D. program.The distribution of the subjects among the caste groups is shown in Table 3. The data36Caste Men WomenBrahmins 75 75Kapu 275 220Balija 153 128Mala 63 99Christians 28 11Muslims 17 51Totals 687 614Table 3: Subject Distribution among Caste Groupscollected included measurements of three groups of variables. Henceforth we shall referto these groups as Group A, Group B, and Group C; Table 4 lists the variables in thesegroupings. Group A consists primarily of pulmonary functions, Group B consists mainlyof anthropometric measurements, and Group C is comprised solely of anthropometricindices. For each subject, the caste to which the subject belonged was also recorded.For each sex within each endogamous group, Dr. Reddy originally proposed to modelseveral of the pulmonary functions, namely vital capacity, forced vital capacity, andforced expiratory volume, using a selection of Group B variables as the explanatoryvariables. This would allow her to determine whether different model parameters werenecessary for some groups. In addition, she intended to model these three pulmonaryfunctions for men and for women, ignoring the groupings, in order to compare her overallmodel with other published results.8.2 The Data Analysis ProblemDue to a set of unusual circumstances, all the original data were lost; however, a largenumber of summary statistics was retained, since the tables containing these statis-tics had been published in the original thesis (see Reddy, 1982). Dr. Reddy consultedSCARL to see if it was possible to reconstruct the data from these statistics, or if thestatistics could somehow be used to do multiple linear regression. If it were possibleto do regression, she was interested in checking versions of the linear model with some37Group A Variables Group B Variables Group C VariablesTidal Volume Age Bi-acromial Breadth IndexVital Capacity Sitting Height Vertex Relative Thoracic IndexInspiratory Capacity Height Vertex Pondoral IndexInspiratory Reserve Volume Bi-acromial Breadth Robusticity IndexExpiratory Capacity Arm Span Pignet Vervaek IndexExpiratory Reserve Volume Chest Circumference (insp.) Nutritional IndexForced Vital Capacity Chest Circumference (exp.) Weight IndexForced Expiratory Volume (FEV 1 ) Body Weight Chest Girth IndexBlood Pressure (diastolic) Body Surface Area Span Stature IndexBlood Pressure (systolic) Chest Width Sitting Height Stature IndexPulse Rate Bi-iliac DiameterBody Temperature Trochanteric DiameterGirth of AbdomenLength of SternumAntero-posterior Chest Diam.Chest VolumeChest ExpansionTrunk Surface Area ITrunk Surface Area IITable 4: Measurements Collectedof the Group C variables (anthropometric indices) as explanatory variables in additionto versions using the Group B variables, since the indices included some importantanthropometric measurements for which separate data were not available.In all, for each sex, she wanted to test six different versions of the basic linear model;these versions would be created by using different combinations of dependent and ex-planatory variables as detailed in Section 11. We shall refer to these different versionsof our basic linear model as "model variants".The following is a list of the available statistics.38For the entire sample:• Range, mean, standard deviation, standard error, and coefficient of variation for:1. all variables grouped by sex2. the first eight group A variables grouped by sex and age group3. all group B variables grouped by sex and age group4. all group C variables grouped by sex and age group• results from t-tests (including t-values) of:1. sex differences for mean values of all group A variables2. sex differences for mean values of all group B variables3. sex differences for mean values of all group C variables4. differences between smokers and non-smokers (men only) 2 for all group Aand group B variables• correlation coefficient estimates for each sex for:1. age versus all group A variables2. age versus all other group B variables3. age versus all group C variables4. the first eight group A variables versus all group B variables5. the first eight group A variables versus all group C variables• inter-correlations of thirteen group B variables for each sex• multiple regression coefficients for vital capacity as a function of 13 group B vari-ables (including age and height) for each sex2 There were no data on female smokers, since they are extremely rare in India.39For each caste group:• Range, mean, standard deviation, standard error, and coefficient of variation forall variables grouped by sexGroup Comparisons:• t-test results (including t-values) by sex for each of the 15 possible pairwise com-parisons between caste groups for:1. mean values of the first eight group A variables2. mean values of all group B variables3. mean values of all group C variablesIn all the tables which showed t-test results, there were some results which were signif-icant at the a-level of .05.To investigate whether a useful model could be checked using the available statistics,we developed the Gaussian distribution based likelihood function which was consistentwith the working assumption of sufficiency of the retained statistics.409 Sufficient StatisticsWe adopt the structural model:yi; = pi + T:ij + fi j= 1,...,6; j = 1,...,ni,where i = group number, j = observation number,xij (1 ) —(q) ^(q)Xij^- X..±.. (Q)and q = covariate number (q =1,...,Q).The stochastic modelling assumptions areeijti H( 0, o-2 ),so thatX =Yij^Ar(pi + xT:ij 13 0-2) •In this model, the value of is the same for all groups; in other words, we are assumingthat each explanatory variable has the same coefficient for all six groups. The interceptterm, pi, on the other hand, can have a different value for each group.This model may seem unrealistic, but since we are only interested in whether one setof parameter values can describe all groups, this model can be used to test pi = pi forall i, j such that i j. If this hypothesis is rejected, we will conclude that the groupscannot be described by the same equation. Note that failure to reject the hypothesiswould not imply that the same equation is appropriate for all groups, since with theavailable statistics, we cannot check to see whether should vary from group to group.41I^/T:il^(xil (1) — ±.. (1) )=X^ (xini (1) —^. . .T:ini9.1 NotationWe shall follow the convention of using uppercase letters for matrices, boldface lowercaseletters for vectors, and unemphasized lowercase letters for scalar quantities.Let n =^ni, and let 1( k) denote a k-dimensional vector of l's. In addition, let usdefine the following vectors and matrices:^= (0 (1)^(Q))^It = (Pi 1-(n/ )^• • , /16 1 (n6 ) )Ei = (Eil , • • • , Eini)E = (E 1 ,..., E6/XT:ij =^(1)— X..(1)^x . . (Q)^±. . (Q))(xil(Q) - .. (Q ) )(xini (Q) — ±..(Q))where Xi is a matrix with dimensions n i x QX =I X6where X is a matrix with dimensions n x QVi = (Vil , • • • , Yin; )= (Y1 / , • • • , y6 / )= — g.. 1 (n)-429.2 Calculation of Sufficient StatisticsUsing the notation defined in Section 9.1, our model can be expressed as:y= A-FX, 3+E.Before finding the sufficient statistics for this model, we shall prove the following lemmaswhich the author has derived.Lemma 16 niEE(yii — — X 7,/,ijO ) 2i=1 j=1ng! — 2y 'it+211 /X13+ 11'11-2VX,3+ O ix'xs.Proof:6 niEE(yii — pi — xTi: iao) 2 =^— — Xi=1 j=1= (y - - x0)'(y - - xi3 )^y ^Wit+ O'X'Xi3— 2y 'it— 2y 'XO-F2it'XI3^= y^11 '11 + 'X'X,Q— 2y 'tt— 2y 'X,3 2y..1 (n' )X0 2/./. 'Xi(3since'(n) X = [0].Thus6^nn sEE (^— /2. —^) 2 =^s'x'xgi=1 j=1— 2y 'it— 2 "Y^+ 211 'X,343=^+^ Wit+ gix'xs- 2y^4 IX/3 + 2p 'X/3since' = E?-1 E;t1(yii Y..) 2 = E?-1^y!.; -^=^ne ■Lemma 26y = EProof:Y^= Y 11 1 1-(ni) + Y2' P2 1-(n2)^• • + Y6' µ6 1 (n6 )= µ1y1 ' 1 (71) + P2Y2 1 (n2) + • • • + µ6y6 ' 1 (76)=^A272h.^fisnds.6=Now we can proceed to find the sufficient statistics for our model. The unknown pa-rameter O = (11,0,o- 2) and1f(Yi.iie) =^1 exP{-T;T(y1.; - - x;,,,O) 2 }.By using the result of Lemma 1, we find that6 ni13^f(y le) = (27r0.2)-71/2explE E^_ pi -^iv)2cr2 i=1 j=1, 1 „= (27ro-2 ) -ni2 exPt^'i/^+ 'tt /3'X'X/3— 2y — 2 -Y 'X,3 +2A'Xi311.44Therefore,f (y I 0)1^1^_2^1^,exp{--721-log 2ircr2 — 20.2 y y — 20.2 ny — 2 u uo-2 "1^1^Q2^1— 2a2 f3 =-1/ — ;VI 'Xi3 }cr' r^1^1^_2^1 v—,6^1y= exPi — Tr.2^—^n^--2- 2_, pinigi.^X02Q2 CT1^ 1^1— —log 27ra2 — -7A , X0 — 207A , A T—(72 0„X X13},where we have substituted the result of Lemma 2 in the last equation.We find from this expression that (i) '""Y; e; yi., i = 1, . . . , 6; 'X) are sufficient statisticsfor the family of distributions indexed by 0.Now we must see whether these statistics can be obtained from what we have at ourdisposal.As noted in Section 8.2, a number of statistics are available, and we wish to check sixdifferent model variants for each sex. For simplicity, assume for now that we wish totest one of the model variants for males, where y is vital capacity, for example, and wehave Q explanatory variables (x's). From the tables at our disposal, we can constructthe following matrices for this model variant. The subscripts, which will be dropped inlater references, indicate the dimensions of each matrix. These matrices are:Q+ 1 ) X( Q+ 1 )Sdiagonal(Q+ 1 )x(Q+ 1 )M6 X(Q-F1)Ndiagonal6x6the correlation matrix for the x's and y for all males;the matrix of standard deviations of the x's and y for all males;the matrix of means of the x's and y for each group;the matrix of the group sizes.45From these matrices, we calculate the sufficient statistics as follows. First, gi., i =1,^, 6 is obtained directly from M. Now, let n = tr(N). To find g.., calculateY±-.. (1)( 146)NM (Q)y..Next, calculate (n — 1)SRS, which can be partitioned thus:[X'X X' f/"'Y 'X -Y iy •This gives us "U and I/ 'X.We have shown here that we do indeed have sufficient statistics for the model fam-ily described at the beginning of Section 9. In the next section, we shall develop anappropriate testing procedure for the substantive problem addressed by this work.4610 Derivation of the Method of AnalysisSince our goal is to compare regression model intercepts, that is, the values of pi, wemust calculate their estimates as well as estimates of their expected values and variances.10.1 NotationLet us define the scalars, vectors and matrices we require, using the conventions estab-lished in Section 9.1 (j3 and are as defined previously in Section 9.1. They areincluded below for completeness.):= (P (1) ,•••, # (')'XT:ij =-7 (X (1) — .. (1) ^(Q) —^(Q)) /°G:i = (x i. (1 ) -^-^/xw :ij^(xij ( 1 )^(1),^xij (Q)^(Q))x i,—— x i (q)Üi.i =^Y.Eij = Eij —^6 ni^ 6 ni^ 6 ni= E E = E E = E E ..xW:yx^ c W:ex^Eli° w:i;^Cw:sx^i=1 j=1 i=1 j=1 i=1 j=1where Cw:sx is a symmetric matrix with dimensions Q x Q. In addition, we define thescalar6 niCW:yy = Ei=1 j=14710.2 Derivation of EstimatorsAs described previously in Section 9, our model is:Yii =^+ z 7!:i.d3 cii ,i =1,...,6; j =1,...,n i ,where i = group number, and j = observation number. We find the least squaresparameter estimators by minimizing the residual sum of squares (RSS) with respect tothe unknown model parameters. This leads to the following estimating equations:aRss^n,^ = —2 E ( yij — pi —^o,aRss^6 n,ago)^ — 2 E E — .0^= o,i=1 j=1where6 niRSS = E E ( yij — pi — x '../3) 2 .T: s9Now, ORSSIap i = 0 entails E .7.L i (yii — pi — xTiiif3) = 0, implying y i . —^xT"ii0 =nipi. From this equation, we obtain(1)(2)It follows that6 n,^ n,RSS = E E [ 6yii — (Yi. — Gx I ,o) x^>: E^— x^)2.i=1 j=1 i=1 j=1We then find that ORSS/0 /3 ( q ) = 0 implies6 n,^ 6 n,^EE^= E E zwl :Jig1.1 j=1 i=1 j=1for every coordinate q. This gives us the equation c^= Cw.t4. Thus11=48i=1 j=1Now, sinceN.; = +^+andx G'.0 + 4 - ,we may subtract Equation (4) from Equation (3), giving yid = xwl :ijo + Eii. Then=cw:^Cw:xx 0+cw:ex . Thus, from Equation (2), ,,j = C-1 c^= C-1 (c w:xxr-A+cw ).ys W:xx w:Yx^w:sx‘^:ex/We conclude that= 0 + C- w:tz^ (5 )Returning to Equation (1), we obtain pi = yi. — xL13. On substituting Equation (4)into Equation (1), we get^=^-I- x 0'. 13^— x G'.i;e1. Equation (5) then gives= + x ' + Ei. — X (0 + C -1 C ). Finally,G:i^G:i^wasW:exfti =^—^C -1^.G:i W:xx •• ' cx(6)^Equation (6) gives us E(fti) =^A i is unbiased. Next, we need to find the variance ofA i . To that end, we find from Equation (6) that — = Ej. — x ' C -1 c . SoG:i^ :xx W:cxVar(ai) =^—Thus,= Var(ei. — x C c )G:i W:xx W:ex= Var(zc )— 2 Cov(x^,G:i W:xz W:cx^G:i W:xx •• :CTVar(•).Var(ii i ) = x '.C -1 Cov(c )C' x — 2 Cov(xL wC- 1.xc w:ex , ei.)^Var(ei.). (7)G:s W:xx^W:ex^G:iTo simplify this expression, we need the following lemmas, derived by Dr. Greig andexpanded by the author.(3)(4)49Lemma 3COV(X 1 C -1 CG:i W : =. W:es Ei- ) = 0.Proof:For simplicity, let x = c^y = '4., and a'= x^. Note that a' is a 1 x Q vector.G:i W:xxThenCov(x^c^,^= Cov(a'x,y)G:i w:xx W:esQ= Cov(Ea qxq ,y)q=1QE cov(aqx q , y)q=1QE aq C ov (x q , y).q=1But6 niCov(x q , y) = Cov[ (E E E-,., x-,(0) ,i=1 j=16 ni= E Ei=1 j=1Now, sinceniCOV(EijEj .) = COV(Eij 1/ni E eik) — Var(EL),k=1we find that, due to independence of the fii's, Cov(Eii, ei.) = ni^ni—^= 0 for all q. Itfollows that EQq_ i a q Cov(xq , y) = O. •50Lemma 4COV(C^= 0-2C^w:cx^w:xz •Proof:Let c( g) denote the et element of cw:es, and C(P9) the entry in the pth row and gillw: x ^W:zzcolumn (and the eh row and pth column) of Cw:xx• Then6 niE(c(q ) ) = E (E EW : Ex i=1 j=16 niE E ivq)E(Eii )i=1 j=1= 0.w:E.Thus Var(c^) = ERc( q) ) 2] and Cov(c ^,c(P) ) = E(c( g) c(P) ). So Cov(c^) =_-w:e.^w:E.^ w:e. w:exw:e. w:e.E (e w:"c w' : ex ). Now,6 niw:E. = E E (6 ij(0)i=1 j=1since^6 74^6E E = E Ei E = 0.^i=1 j=1^1^i=1 j=1Thus6 ni^6 nkE(c( g) c(1') )W:ex W:es = E [ (E E ii.i(q)) • (E E coki(P))1=1 j=1^k=1 1=1^6 ni^ ni nk^= E [ E E (fi.0 2 .i.;(q)i.j(P)J+E[EEE^iii(qki(P)j=1^ i0k j=1 1=16+ E [E E (Eij ci,)i=1 j016 ni= E E1=1 j=1since expectation of the last two terms is zero due to independence of the cii 's. Itfollows that E(c( q) c(P) ) = o-2 C (Pq ) , and in turn that E(c^c ' ) = cr2Cw:ex w:cz^Thusw:ex^w:x.^ w:x.•Cov(c,„) = o-2Cw:rr. ■51Now, substituting the results of Lemmas 3 and 4 into Equation (7), we geta2Var(iii)^x c-i cr2cG:i W:xx^W:xx C -1 X .^w,xz G:.^ni22^C-1x:zz G:= a  G:i wsAlso,^Cov(iii, [Li) = Cov(fzi —^— pj)^= C OV[(ei. — X GIS w- 1x.Cw:f. )^— X GIiC= C 00[(X C —1 C ) (X C —1 C )]G:i w :xx W:ex^G:j w :sx W:ex= E(x C -1 c x C -1 cG:i W:xx W:ex G:j W :xx W:ex)^= E(x C -1^I C 1 X )G:i w :xx Wsex W:ex w :zz G:3= X C -1 C 00(C )C1G:i^W:ex, vv :zs X 0:3= x C -1 0-2C C-1 x^G:i W:xx^W:xx W:xx G:j2 I C-1.= Cr XG:i w:zr G:jThenVar(fli ji) = Var(iii) Var(fL) — 2 Cov(it i , Ai )2^ u2Cr^= (0.2x i .c-i x + __) (0.2x .C-1 x^__) _ 2 (0.2x r .c-iG:s W:sx G:i^ W:xx 0:3 G:i W:xx G:j)ni ni^a2 (x C' x + — + x 31 ^x^—1 — 2x C -1 x )G:i W:sz G:s^ w:x. 0:i G:i W:xs G:3n i^nj1^1=G:i^W:xx^•-••1^G.352Since o2 is unknown, we must estimate it. Let f denote the number of degrees offreedom in our model variant. This gives us f = 6 4- Q. Then2CI1^6 niE E(Yii — Pii)21^6 nin f E^—^— 0Z:1^XT:ij3) —^13 } 22=1 3 =11^6 nin^f^1(gij^41,:ij44)2j=6 ni1[ EM.02n — f 6 ni—2 E E1=1 j=1^36 ni+ E E s-1 c„,^c-1 cw:xx ..:y. W:xx W:yx] •1=1 j=1The last equation is obtained by substituting Equation (2). Therefore,2^1w: yx W:xxCW:Yx(^ -1— 2 c Cn — f cw:YY6 ni+ EEC C- 1^..x C-1W:Yr W:xx W33 W:i1 W:xx W:Yx )1=1 j=1n1f (c—^W•yy— 2c C -1 cw:yx W:xx w:y.c^C C-1 C,„W:YX W:XX^:XX W:XXsince EL.. 1 E'I=i^= Cw:... Finally,crA 2 1n — f (CW:yy — c^cW:Yx W:sx W:Yx ). (8)5310.2.1 A Restricted ModelIn this section we consider the model of Section 9 with the additional restriction that= p for all i:yi; = +^+To distinguish the estimators from the two models, our stochastic modelling assumptionsin this case areUsing the notation defined in Section 9.1 and proceeding in the manner of the derivationof Equations 1, 2, and 8 in Section 10.2, we find that= g..igo = (x xT:s,^T:g3andr / (70 = [Y Y Y^—1 xT:iiY)/(n Q - 1). (9 )We shall refer to this model in Section 10.4 when we are considering hypothesis tests.10.3 Calculation of EstimatesReturning to our less restricted model, and using 172 in place of u2 , we find that theestimators needed for our analysis are:gi.^=^c,„„. .For these estimators,2V ar(fti) = 6-2x C -1 x . —niG:i w:zz G:tand2 / ry-1= A XG:i‘-'w:xxXG:i'where ef 2 =^— c C -1 c ).n —f W:YY^W:ys W:xx W:Yx54To find the estimates, we need to calculate the following quantities: yi.; ni; 7 C^•W:xx 7c^• and c^.W:yx 7^W:yyAs in Section 9.2, we seek to test one of the model variants for males, and construct thefollowing matrices for this variant. We repeat the matrices here for convenience:AQ+ 1 )X(Q+ 1 )Sdiagonal(Q+1)X(Q+1)Al6x(Q+1)Ndiagonal6 x6the correlation matrix for x's and y for all males;the matrix of standard deviations of x's and y for all males;the matrix of means of x's and y for each group;the matrix of group sizes.From these matrices, we can calculate the statistics we require. First, yi. and ni areobtained directly from M and N respectively. Now, let n = tr(N). To find x0 ,,,calculate the augmented matrixG6x(Q+1) M itaLLVInXG:1XG:21XG:3XG:4XG:5XG:6 Next, let CT = (n-1)SRS = the total sum of squares matrix, and let CB = G'NG = thematrix of between group sums of squares and cross-products. Then Cw = CT — CBthe matrix of within group sums of squares and cross-products. Note that Cw is a(Q + 1) x (Q + 1) matrix.It can be partitioned thus:[ Cw:rr CW:yzThis gives us Cw., cw,, and c yy •CW:yx CW:yy55Note that, for computational convenience, in the actual construction of the matricesfor each sex, we included all the x's and y's which were to be considered in any of thesix model variants for that sex. Then, when calculating the estimates for a particularvariant, we simply selected the appropriate rows and columns from the matrices.10.4 Testing MethodologyOnce we have calculated the necessary estimates for each of the twelve model variantswe propose to check, we wish to test (for each variant) whether the intercepts of theregression planes for all six groups are the same, i.e. our null hypothesis is:Ho^= /12— =To test this hypothesis we could perform the F-test, which some would consider theclassic test. To implement it, we letS — (n — Q^(n — Q — 6)er25&2where 6- 2 and c-rg are given by Equations 8 and 9 respectively. Under the assumptionsof the restricted model described in Section 10.2.1, the statistic S has an F5,n_Q_6distribution and may be used as a test of the null hypothesis Ho (see Graybill).However, in the event that there was not equality of pi 's, Dr. Reddy wished to be ableto say which pairs of pi's differed. We therefore propose to test 1/ 0 by using fifteenpairwise tests, and adopting the hypothesis^Hof : pt = j i = 1,...,6; j = 1,...,6; i^j.Since under the assumption of normality,—NIV ar(P i —we utilize this as our test statistic. We calculate this statistic for each of the pairwisetests, and obtain the p-values for these statistics using the t distribution. Let p i , ... , p15be the fifteen p-values which result from these tests as shown in Table 5.has a Student's t distribution,56Groups Compared p-value1 versus 2 P11 versus 3 P21 versus 4 p31 versus 5 1341 versus 6 P52 versus 3 P62 versus 4 P72 versus 5 P82 versus 6 P93 versus 4 Plo3 versus 5 P113 versus 6 P124 versus 5 P134 versus 6 P145 versus 6 P15Table 5: Pairwise Comparisons for Testing Equality of ,ui 's57If we were simply testing to see whether two particular groups have the same intercept,with an a-level of .05, we would just check to see if the p-value resulting from thatparticular pairwise comparison were less than .05. If it were, we would reject the nullhypothesis of equality. Note that this would be equivalent to checking only one of thefifteen paired tests included in N. In this case, however, we decided, in consultationwith Dr. Reddy, to test H!, (rather than H0 ) with an a-level of 0.05.We can do this by using the modified Bonferroni test known as Holm's sequentiallyrejective Bonferroni procedure. A brief development of Holm's sequentially rejectiveBonferroni procedure follows; for more complete details see Hochberg and Tamhane(1987).The standard Bonferroni test is based on the Bonferroni inequality, which states thatP(R1 U R2 U . . . R,,) <^1 P(R1). The standard Bonferroni test guarantees the leftside is < a by making P(Rj) a/n for every i. If we were to use this test in thecurrent situation, we would compare each of the p-values from the pairwise comparisonswith .05/15. If any of these p-values were less than .05/15, we would reject the nullhypothesis, Ho , with an a-level of .05. Note that the standard Bonferroni test is asingle-step procedure, and requires that all test results be checked.We would prefer instead to use a step-down procedure, since it will be more powerful.In general, a step-down procedure begins by testing the overall intersection hypothesisand then steps down through the hierarchy of implied hypotheses. If any hypothesis isretained, all of the implied hypotheses will be retained without further tests. In thisway, a hypothesis is tested if and only if all of its implying hypotheses are rejected.A general method for step-down procedures, called the closure method, was proposed byMarcus, Peritz, and Gabriel (1976). For closed testing procedures, we let H1(1 < i < 7n)be a finite family of hypotheses, and we form the closure of this family by taking allnonempty intersections Hp = niEpHi for P C 1, 2, ... , m. If an a-level test of eachhypothesis Hp is available, then the closed testing procedure rejects any hypothesis Hpif and only if every HQ is rejected by its associated a-level test for all Q D P.Note that the number of tests in a closed testing procedure increases exponentially withm. It is therefore desirable to develop a shortcut method. One such shortcut is Holm's58sequentially rejective Bonferroni procedure, which is carried out in our case by followingthese steps:1. We arrange Pi,^, P15 in non-decreasing order p( i ), • • • ,p(15).2. At the same time, we record the order (in size) of /3 1 ,^,Pis . From this operation,we get a list of values k, k = 1,... ,15 such that pk = p(j) for some j, j^1,^, 15,so that later we can identify the pairing which produced each p(j).3. We calculateqi ,^, qi5, where qi =^, q2 = it, • • • ,q1.5 =4. Beginning with p( i), *we compare p(j) with qj. If p(j) is less than qi , we check thevalue of k to find the associated pair, reject H,C for this pair, increment j by one,and repeat from *. If p(j) is greater than qj, we store the current value of j in1, discontinue testing, and retain the null hypothesis 14 for all pairs associatedwith p(j) where j > 1. This can be done conveniently by examining the signs of(qj — pm), j^1,^, 15. If the sign is positive, p(j) is less than qj.5. At the end of this procedure we can say, with an a-level of .05, that all of the pairswith an associated p(j) where j < 1 are significantly different.The process of calculating the test statistics, conducting the pairwise tests, and perform-ing the modified Bonferroni test as detailed above is performed for each of the twelvemodel variants.In future work, we might consider using more recent developments in the theory of mul-tiple comparisons. However, given our uncertainty about the validity of the underlyingworking assumptions pending further exploration by the investigator, we decided not topursue further refinement of the analysis at this time.5911 ResultsThe variables of interest are shown in Table 6. The caste groups are:Group 1: BrahminsGroup 2: KapuGroup 3: BalijaGroup 4: MalaGroup 5: ChristiansGroup 6: Muslims.The basic model as specified in Section 10 was tested for twelve different model variants,six for men and six for women. These model variants are listed in Table 7. For modelvariants m.4, m.5, m.6, w.1, w.4, w.5, and w.6, there appeared to be no significantdifferences between the intercepts for the six groups. The results from the modifiedBonferroni test for the other model variants are summarized in Table 8 and depictedgraphically in Figure 11.Tables 9 and 10 display the results of both the individual pairwise tests and the mod-ified Bonferroni test for all of the model variants which showed any significant groupdifferences. The "I" columns contain an asterisk if the indicated pairwise comparisonresulted in a p-value of less than .05. We can say that each of these pairs is significantlydifferent with an a-level of .05. The "B" columns contain asterisks for those pairwisecomparisons with an associated p(;) where j < 1 (where those terms are as defined inSection 10.4). For the indicated model variant, we can say that, taken as a group, allthese pairs are significantly different with an a-level of .05.Appendix I contains the Splus computer routines which were used for the more complexcalculations.For all twelve model variants tested, the actual p-values, p3 (as defined in Section 10.4)can be found in Appendix II. Table 13 in the same appendix shows the detailed resultsof the modified Bonferroni test.60Variable I DescriptionGroup B VariablesX1X2X3X4X5X6X7AgeHeight VertexSitting Height VertexArm SpanBody WeightChest ExpansionBi-acromial BreadthGroup C VariablesY1Y2Y3Y4Y5Y6Y7Y8Y9Y10YllY12Y13Body Surface AreaTrunk Surface AreaChest VolumeBi-Acromial Breadth IndexRelative Thoracic IndexPondoral IndexRobusticity IndexPignet Vervack IndexNutritional IndexWeight IndexChest Girth IndexSpan Stature IndexSitting Height Stature IndexPulmonary FunctionsZ1Z2Z3Vital CapacityForced Vital CapacityForced Expiratory VolumeTable 6: Variables of Interest61VariantNumberDependentVariableExplanatoryVariablesMenm.1m.2m.3m.4m.5m.6ZiZ2Z3ZiZ2Z3X1 to X5, X7X1 to X5, X7X1 to X5, X7X6, Yi to Y7, Y11 to Y13X6, Yi to Y7, Yii to Y13X6, 1/1 to Y6, Y8, Yll to Y13Womenw.1w.2w.3w.4w.5w.6ZiZ2Z3ZiZ2Z3X1 to X5, X7X1 to X5, X7X1 to X5, X7X6, 1/1 to Y7, Yll to Y13X6, Y1 to Y7, Yll to Y13X6, Y1 to Y7, Y11 to Y13Table 7: Model Variants Tested62Groups ComparedResultsVariant m.1 Variant m.2 Variant m.3 II Variant w.2 Variant w.31 versus 2 * * *1 versus 3 * * * *1 versus 4 * * * *1 versus 5 *1 versus 6 * * * * *2 versus 3 * * * *2 versus 4 * *2 versus 5 * *2 versus 6 * * * *3 versus 4 * * *3 versus 5 * *3 versus 6 * * * *4 versus 5 *4 versus 6 * * * *5 versus 6 * * *Table 8: Test Results from Bonferroni Comparisons for Men and Women"*" indicates the pairs in each variant which are significantly differentwith an a-level of .05 for the entire model variant.63654• 9• 3Model Variant m.11 Group Key 1. Brahmins2. Kapu3. Balija4. Mala5. Christians6. MuslimsModel Variant m.21Model Variant w.216365Model Variant m.31Model Variant w.31Figure 11: Results for Model Variants with Significant Differences at a-level of .05.Points not connected to each otherrepresent groups which are significantly different at this level.64Groups ComparedVariant m.1 Variant m.2 Variant m.3I B I B I B1 versus 2 * * * *1 versus 3 * * * * * *1 versus 4 * * * *1 versus 5 * *1 versus 6 * * * * * *2 versus 3 * * * * * *2 versus 4 * *2 versus 5 * * * * *2 versus 6 * * * * * *3 versus 4 * * * * * *3 versus 5 * * * *3 versus 6 * * * * * *4 versus 5 * *4 versus 6 * * * * * *5 versus 6 * * * * * *Table 9: Test Results from Individual and Bonferroni Comparisons for Men"I" indicates individual pairwise comparisons."B" indicates modified Bonferroni test for the entire model variant."*" indicates significance at the level of a = .05.65Groups ComparedVariant w.1 Variant w.2 Variant w.3I^B I B I1 versus 21 versus 31 versus 4** *******1 versus 51 versus 62 versus 3* ******2 versus 42 versus 52 versus 6* ****3 versus 43 versus 53 versus 6** *4 versus 54 versus 65 versus 6* *Table 10: Test Results from Individual and Bonferroni Comparisons for Women"I" indicates individual pairwise comparisons."B" indicates modified Bonferroni test for the entire model variant."*" indicates significance at the level of a = .05.6612 ConclusionsIt is important to note that the absence of the original data forced us to make a numberof unverifiable working assumptions when we were building our model. We assumed thesuitability of a linear model with a Gaussian-distributed error term with equal variancefor all observations. In addition, we assumed a common value of /3 for all groups. Note,however, that the latter assumption, even if violated, does not invalidate our procedure.The assumption simply renders it conservative, as discussed in Section 9.Data quality is always a concern, but it is doubly so in our case. Not only might there bemisgivings about the method of subject selection, but in addition, the accuracy of thesummary statistics themselves could be questioned since errors might have been madein their production and transcription. In the absence of the data, the reliability of thesetables cannot be checked.Nonetheless, the methodology which was developed to handle this problem appearssound, and it is remarkable for its simplicity and its ability to produce an analysis froma minimal amount of information.67References for Part IIGraybill, F.A. (1976). Theory and Application of the Linear Model. North Scituate,Massachusetts: Duxbury Press, 189-190.Hochberg, Y., and Tamhane, A.C. (1987). Multiple Comparison Procedures. New York:Wiley, 53-58.Marcus, R., Peritz, E., and Gabriel, K.R. (1976). On closed testing procedures withspecial reference to ordered analysis of variance. Biometrika, 63, 655-660.Reddy, L. (1982). Study of Pulmonary Functions and their Relation to AnthropometricParameters in Men and Women of Chittoor District of Andra Pradesh. Ph.D.Dissertation, Department of Anthropology, University of Madras.68Appendix I: Splus Routines for Part IIThe following are the computer routines (written in Splus) which were used for thevarious calculations.Calcint.s: used to calculate ji and Cov(#) for each of the six model variants for*sd.m%*%corr.m%*%sd.mJ16matrix(1,1,6)ybar_a/n)*J16%*%size.m%*%mean.mJ66_matrix(1,6,6)zl mean.m - (1/n)*J66%*%size.m%*%mean.mwlt1 - t(zl)%*%size.m%*%zlwg_wl[indvector,indvector]w_wl[indvector,i]omegawl[i,i]zg_z1[,indvector]zz1[,i]betahat_solve(wg)%*%wsigma2hat_(1/(n-f))*(omega - t(w)%*%solve(wg)%*%w)grandmean(1/n)*J16%*%size.m%*%mean.m[,i]ybarbarrep(grandmean,6)muhat_ybarbar + z - zg%*%betahatsigmamuhat_sigma2hat[1,1]*(solve(size.m) + zg%*%solve(wg)%*%t(zg))n =i =indvector =f =total number of males studied.the number of the column (and row)of the standard deviation and correlation matrices whichcontains the dependent variable for the current model variant.the vector of column (and row) numbersfor the explanatory variables for the current model variant.the number of degrees of freedom for the current model variant.For each model variant, the above values are set, the routine is called by the"source" statement, and then the values of "muhat" and "sigmamuhat" are storedin the appropriate storage registers for the current model variant. The followingexample shows the sequence for variant m.1.n 687i —21iiidvector_c(1,2,3,4,5,7)f_ 12source("calcint.s")muhat.m1 muhatsigmamuhat.ml_sigmamuhat69Calcint.sw: used to calculate # and Cov(#) for each of the six model variants*sd.w%*%corr.w%*%sd.wJ16matrix(1,1,6)ybar_a/n)*J16%*%size.w%*%mean.wJ66_matrix(1,6,6)zl mean.w - (1/n)*J66%*%size.w%*%mean.wwlt1 - t(z1)%*%size.w%*%z1wgwl[indvector,indvector]w_wl[indvector,i]omega_wl[1,1]zg_z1[,indvector]zz1[,i]betahat_solve(wg)%*%wsigma2hat_(1/(n-f))*(omega - t(w)%*%solve(wg)%*%w)grandmean(1/n)*J16%*%size.w%*%mean.w[,i)ybarbar_rep(grandmean,6)muhat_ybarbar + z - zg%*%betahatsigmamuhat_sigma2hat[1,1]*(solve(size.w) + zg%*%solve(wg)%*%t(zg))n = total number of females studied.i = the number of the column (and row)of the standard deviation and correlation matrices whichcontains the dependent variable for the current model variant.indvector^the vector of column (and row) numbersfor the explanatory variables for the current model variant.f = the number of degrees of freedom for the current model variant.For each model variant, the above values are set, the routine is called by the"source" statement, and then the values of "muhat" and "sigmamuhat" are storedin the appropriate storage registers for the current model variant. The followingexample shows the sequence for variant w.l.n6141 21indvector_c(1,2,3,4,5,7)f 12source("calcint.sw")muhat.wl_muhatsigmamuhat.wl_sigmamuhatmuhat.wl70testloop.s: used to calculate the test statistics and probability values for each of thefifteen pairwise tests for each model variant.testvectorrep(0,15)test index 1for (i in 2:6)testvector[testindex]_(mu[1]-mu[i])/sqrt(sigma[1,1]+sigma[i,i] - 2*sigma[1,i])testindex_testindex+11for (i in 3:6)1testvector[testindex]_(mu[2]-mu[i])/sqrt(sigma[2,2]+ sigma[i,i] - 2*sigma[2,i])testindex_testindex+1for (i in 4:6)testvector[testindex]_(mu[3]-mu[i])/sqrt(sigma[3,3]+ sigma[i,i] - 2*sigma[3,i])testindex_testindex+11for (i in 5:6){testvector[testindex]_(mu[4]-mu[i])/sqrt(sigma[4,4]+ sigma[i,i] - 2*sigma[4 i i])testindex_testindex+1}testvector[testindex](mu[5]-mu[6])/sqrt(sigma[5,5]+ sigma[6,6] - 2*sigma[5,6])pvector_2*(1-pnorm(abs(testvector)))mu = # for the current model variant.sigma = Cov(it) for the current model variant.testvector = the 15 x 1 vector of test statisticsfor group comparisons in the order indicated in Table 5.pvector = the 15 x 1 vector of p-valuesfor these test statistics.71Appendix II: Test Results for Part IIGroups Compared p-valueModel Variantm.1 m.2 m.3 m.4 m.5 m.61 versus 2 Pi. 0.00 0.12 0.00 0.62 0.09 0.341 versus 3 p2 0.00 0.00 0.00 0.34 0.06 0.361 versus 4 p3 0.01 0.26 0.00 0.07 0.08 0.941 versus 5 p4 0.36 0.35 0.02 0.96 0.51 0.301 versus 6 p5 0.00 0.00 0.00 0.13 0.08 0.122 versus 3 P6 0.00 0.00 0.00 0.51 0.70 0.982 versus 4 p7 0.72 0.90 0.00 0.08 0.64 0.402 versus 5 Ps 0.00 0.04 0.00 0.78 0.66 0.612 versus 6 p9 0.00 0.00 0.00 0.18 0.33 0.263 versus 4 Pio 0.00 0.00 0.00 0.23 0.86 0.433 versus 5 Pit 0.00 0.06 0.00 0.54 0.54 0.643 versus 6 P12 0.00 0.00 0.00 0.29 0.42 0.274 versus 5 P13 0.01 0.08 0.24 0.18 0.50 0.344 versus 6 P14 0.00 0.00 0.00 0.75 0.51 0.145 versus 6 P15 0.00 0.00 0.00 0.20 0.28 0.55Table 11: p-values from Pairwise Comparisons for Testing Equality of p i 's for Men72Groups Compared p-valueModel Variantw.1 w.2 w.3 w.4 w.5 w.61 versus 2 Th. 0.04 0.31 0.00 0.05 0.12 0.181 versus 3 P2 0.29 0.09 0.00 0.10 0.13 0.121 versus 4 p3 0.80 0.00 0.00 0.71 0.76 0.911 versus 5 P4 0.75 0.77 0.02 0.55 0.55 0.051 versus 6 p5 0.19 0.00 0.00 0.16 0.50 0.132 versus 3 p6 0.26 0.31 0.00 0.74 0.99 0.742 versus 4 1)7 0.01 0.00 0.07 0.08 0.17 0.172 versus 5 p8 0.56 0.89 0.98 0.79 0.93 0.142 versus 6 p9 0.81 0.00 0.36 0.90 0.54 0.553 versus 4 Pio 0.15 0.03 0.35 0.18 0.20 0.123 versus 5 p" 0.86 0.63 0.25 0.88 0.93 0.183 versus 6 P12 0.60 0.00 0.20 0.91 0.56 0.754 versus 5 P13 0.66 0.16 0.46 0.67 0.65 0.054 versus 6 P14 0.11 0.00 0.63 0.25 0.66 0.135 versus 6 P15 0.67 0.05 0.65 0.85 0.84 0.28Table 12: p-values from Pairwise Comparisons for Testing Equality of it i 's for Women73Variant m.1 Variant m.2 Variant m.3 Variant w.2 Variant w.3P(i)^ii Sgn k Sgn k Sgn k^fi Sgn k Sgn kP(1) + 5 + 12 + 1 + 14 + 1P(2) + 6 + 15 + 2 + 12 + 3P(3) + 9 + 9 + 5 + 9 + 5P(4) + 12 + 5 + 6 + 3 + 6P(5 ) + 14 + 14 + 9 + 7 + 2P(s) + 15 + 6 + 10 + 5 - 4P(7) + 10 + 10 + 11 - 10 - 7P(8) + 2 + 2 + 12 - 15 - 12P(s) + 11 - 8 + 14 - 2 - 11P(lo) + 8 - 11 + 15 - 13 - 10P(n) + 1 - 13 + 7 - 1 - 9P(12) + 13 - 1 + 8 - 6 - 13P(13) + 3 - 3 + 3 - 11 - 14P(14) - 4 - 4 + 4 - 4 - 15P(15) - 7 - 7 - 13 - 8 - 8Table 13: Test Results from Bonferroni MethodA "+" in the "Sgn" column indicates j < 1 for the associated P(j),where p(i) and 1 are as defined in Section 10.4.For each variant, all pairs (considered as a group) with a "+" in this columnare significantly different at the level of a = .05."k" is the value for finding the associated pk (and pairing) from Table 5.74


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items