A mean squared error of prediction approach to the analysis of the combined array

MISSING IMAGE

Material Information

Title:
A mean squared error of prediction approach to the analysis of the combined array
Physical Description:
ix, 149 leaves : ill. ; 29 cm.
Language:
English
Creator:
O'Donnell, Eileen M., 1967-
Publication Date:

Subjects

Genre:
bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )

Notes

Thesis:
Thesis (Ph. D.)--University of Florida, 1994.
Bibliography:
Includes bibliographical references (leaves 147-148).
Statement of Responsibility:
by Eileen M. O'Donnell.
General Note:
Typescript.
General Note:
Vita.

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 002019631
notis - AKK7080
oclc - 32801350
System ID:
AA00003634:00001

Full Text








A MEAN SQUARED ERROR OF PREDICTION APPROACH TO THE
ANALYSIS OF THE COMBINED ARRAY














By

EILEEN M. O'DONNELL


A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY

UNIVERSITY OF FLORIDA


1994











ACKNOWLEDGEMENTS


I would like to express my sincere gratitude to Dr. Geoff Vining for being my ad-

visor and for all the time and effort he put into helping me complete this dissertation.

His patience, encouragement, and guidance during the sometimes slow and frustrat-

ing periods of independent study enabled me to get to the finish line. Dr. Vining has

taught me much more than I could have learned on my own, and I look forward to

continued collaboration with him during the years to come. I would also like to thank

Dr. John Cornell, Dr. Ronald Randles, and Dr. Jack Elzinga for their support while

serving on my committee. Thanks also go to Dr. Brett Presnell for attending my oral

exams. There are many other teachers to whom I wish to express my gratitude for

being role models and mentors during my education, from elementary school through

graduate school.

In addition, I would like to thank my family: Mom, Dad, Steve, Ed, and Tracey,

who are responsible for the person I am today. My parents placed great value on our

education and we have each achieved great success because of it. Special thanks are

also owed to my classmate, Kannan, who at any time would drop whatever he was

doing to help me with my latest "urgent" problem. There are many other friends

who made my years in graduate school memorable. I would also like to thank Carol

Rozear, who has been extremely helpful in dealing with the many nonacademic details

during my years here and for letting me keep my jet skis at her home. Lastly, I wish

to thank David Wolfe for his kindness, patience, and love and for providing me with

much-needed distractions. David will forever be my best friend.












TABLE OF CONTENTS




ACKNOW LEDGEM ENTS ..................................................... ii

LIST OF FIGURES ........................... .... ......................... v

A B ST R A CT ................................................................... viii

CHAPTERS

1 INTRODUCTION .................................................... 1


2 REVIEW OF LITERATURE ......................................... 6

2.1 A Review of Taguchi's Parameter Design Methods ............ 6
2.2 The Combined Array: An Alternative to Taguchi's Crossed
Array .......................................................... 13
2.3 Run Savings and Estimation Flexibility of the Combined Array 17
2.4 Analysis of a Combined Array Experiment ..................... 19
2.5 The Mean Squared Error Criterion ............................. 21

3 ANALYSIS OF THE COMBINED ARRAY
WHEN THE ASSUMED MODEL IS CORRECT ................. 24

3.1 Introduction ..................................................... 24
3.2 Properties of the Estimated Process Variance .................. 25
3.3 Expressions for the Bias and Variance of the Estimated Process
Variance ....................................................... 28
3.4 Comparing Designs Based on the Mean Squared Error Criterion 30
3.5 Results and Discussion .......................................... 40
3.6 Other Measures of Design Performance ....................... 46
3.7 Practical Application ........................................... 53
3.8 Conclusions ..................................................... 55

4 ANALYSIS OF THE COMBINED ARRAY
WHEN THE ASSUMED MODEL IS UNDERSPECIFIED ....... 59

4.1 Introduction ..................................................... 59






4.2 Properties of the Estimated Process Variance .................. 60
4.3 Expressions for the Bias and Variance of the Estimated Process
V ariance ....................................................... 66
4.4 The IMSE Criterion and Other Measures of Design Performance 67

5 SPECIAL CASES OF MODEL UNDERSPECIFICATION ......... 69

5.1 Introduction ..................................................... 69
5.2 Reduced Expressions for the Bias of the Estimated Loss ....... 70
5.3 Exam ple ......................................................... 73
5.4 Results and Discussion ....................................... 75
5.5 Summary and Conclusions ..................................... 113
5.6 Practical Application ................ .......................... 116

6 CONCLUSIONS AND FUTURE RESEARCH ...................... 126


APPENDICES

A INVARIANCE OF MSE ............................................. 129


B BIAS DERIVATION .................................................. 132


C VARIANCE DERIVATION .......................................... 136


D DERIVATION OF COVARIANCE MATRIX ........................ 138


E SAS PROGRAM ..................................................... 143


REFERENCES ................................................................. 147

BIOGRAPHICAL SKETCH .................................................. 149












LIST OF FIGURES


Figure

2.1 Illustration of a Crossed Array ...............................

3.1 16-run Resolution V Fractional Factorial Design .................

3.2 20-run Plackett-Burman Design ..............................

3.3 24-run Plackett Burman Design ..............................

3.4 24-Run Orthogonal Array ...................................


3.5 The Effect of Correlation

3.6 The Effect of Correlation

3.7 The Effect of Correlation

3.8 The Effect of Correlation

3.9 The Effect of Correlation

3.10 The Effect of Correlation

3.11 The Effect of Correlation

3.12 The Effect of Correlation
(Unweighted Weighted)

3.13 The Effect of Correlation
(Unweighted Weighted)


P


the Average Bias (Unweighted Weighted)

the Average MSE ..................

the Weighted Average MSE ..........

the Maximum-Average MSE ..........

the Weighted Maximum-Average MSE .

the Maximum-Maximum MSE ........

the Weighted Maximum-Maximum MSE

the Average Bias and Average MSE
.............................. the Average Bias and Average MSE

the Average Bias and Average MSE


Case 1 (C2) The Effect of Correlation on the Average Bias ........

Case 1 (C2) The Effect of Correlation on the Weighted Average Bias

Case 1 (C2) The Effect of Correlation on the Average MSE .......

Case 1 (C2) The Effect of Correlation on the Weighted Average MSE

Case 2 (N2) The Effect of Correlation on the Average Bias........


age

9

34

38

39

41

43

44

45

48

49

51

52


56


57

77

78

79

80

82






5.6 Case 2 (N2) The Effect of Correlation on the Weighted Average Bias

5.7 Case 2 (N2) The Effect of Correlation on the Average MSE .......

5.8 Case 2 (N2) The Effect of Correlation on the Weighted Average MSE

5.9 Case 3 (C2 x N) The Effect of Correlation on the Average Bias ....


3 (C2 x N)


3 (C2 x N)

3 (C2 x N)


4 (C2 x N2)

4 (C x N2)


4 (C x N2)

4 (C x N2)


- The Effect of Correlation on the Weighted Average
S. . . . . 88

- The Effect of Correlation on the Average MSE ... 89

- The Effect of Correlation on the Weighted Average
S. . . . 90

- The Effect of Correlation on the Average Bias .... 92

- The Effect of Correlation on the Weighted Average
. . . . . 93

- The Effect of Correlation on the Average MSE 94

- The Effect of Correlation on the Weighted Average
. .. .. .. .. .. .. 95


5.17 Case 5 (C2, N2) The Effect of Correlation on the Average Bias .....

5.18 Case 5 (C2, N2) The Effect of Correlation on the Weighted Average
Bias ..................................................

5.19 Case 5 (C', N2) The Effect of Correlation on the Average MSE .

5.20 Case 5 (C2, N2) The Effect of Correlation on the Weighted Average
MSE ...................................................

5.21 Case 6 (C2, C2 x N) The Effect of Correlation on the Average Bias .

5.22 Case 6 (C2, C2 x N) The Effect of Correlation on the Weighted
A average Bias ........................................... ..

5.23 Case 6 (C2, C2 x N) The Effect of Correlation on the Average MSE

5.24 Case 6 (C2, C2 x N) The Effect of Correlation on the Weighted
A average M SE ............................................

5.25 Case 7 (N2, C x N2) The Effect of Correlation on the Average Bias .

5.26 Case 7 (N2, C x N2) The Effect of Correlation on the Weighted
A average B ias. ............................................

5.27 Case 7 (N2, C x N2) The Effect of Correlation on the Average MSE


5.10 Case
Bias

5.11 Case

5.12 Case
MSE

5.13 Case

5.14 Case
Bias

5.15 Case

5.16 Case
MSE






5.28 Case 7 (N2, C x N2) The Effect of Correlation on the Weighted
Average M SE ............................................ 110
5.29 Case 8 (C2, N2, C2 x N, C x N2) The Effect of Correlation on the
Average Bias ......................................... 111

5.30 Case 8 (C2, N2,2 x N,C x N2) The Effect of Correlation on the
W eighted Average Bias ..................................... 112
5.31 Case 8 (C2,N2, 2 x N,C x N2) The Effect of Correlation on the
Average M SE ............................................ 114
5.32 Case 8 (C2, N2, C2 x N,C x N2) The Effect of Correlation on the
W eighted Average M SE .................................... 115
5.33 Case 1 (C2) The Effect of Correlation on the Average Bias and
Average MSE (Unweighted Weighted) ........................ 118
5.34 Case 2 (N2) The Effect of Correlation on the Average Bias and
Average MSE (Unweighted Weighted) ........................ 119
5.35 Case 3 (C2 x N) The Effect of Correlation on the Average Bias and
Average MSE (Unweighted Weighted) ........................ 120
5.36 Case 4 (C x N2) The Effect of Correlation on the Average Bias and
Average MSE (Unweighted Weighted) ........................ 121

5.37 Case 5 (C2, N2) The Effect of Correlation on the Average Bias and
Average MSE (Unweighted Weighted) ........................ 122

5.38 Case 6 (C2, C2 x N) The Effect of Correlation on the Average Bias
and Average MSE (Unweighted Weighted) .................... 123
5.39 Case 7 (N2, C x N2) The Effect of Correlation on the Average Bias
and Average MSE (Unweighted Weighted) .................... 124

5.40 Case 8 (C2, N2, C2 x N,C x N2) The Effect of Correlation on the
Average Bias and Average NISE (Unweighted Weighted) .......... 125











Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment
of the Requirements for the Degree of
Doctor of Philosophy

A MEAN SQUARED ERROR OF PREDICTION APPROACH TO THE
ANALYSIS OF THE COMBINED ARRAY

By

Eileen M. O'Donnell

August 1994

Chairman: G. Geoffrey Vining
Major Department: Statistics

The combined array is an experimental design used for quality improvement within

a robust parameter design framework. Robust parameter design is concerned with

how to reduce the variation of a product's performance. In this framework, the

factors that influence the quality characteristic of interest are of two types: control

factors and noise factors. Control factors are those that are under the direct control

of the experimenter, both in the experiment and in the process. Noise factors can

be controlled for the purposes of an experiment, but are difficult to control in the

process. The goal of a robust parameter design experiment is to select the values of

the control factors that make the product insensitive to changes in the noise factors.

The combined array approach represents an effort to integrate Taguchi's parameter

design principles with more well-established and efficient statistical techniques. This

approach uses a single design matrix and models the response directly as a function

of the control and noise factors. The fitted model is used to determine settings of

the control factors that will minimize an appropriate loss function involving the noise

factors. Often, the most appropriate loss function is simply the resulting process






variance, recognizing that the noise factors, although fixed in the experiment are

actually random effects in the process. When the focus of an experiment is to find

control factor settings that achieve a target with minimum variability, it is vital

to understand the properties of the process variance estimator. This dissertation

represents the first formal attempt at deriving such properties.

The derivation of the properties of the process variance estimator is dependent

upon the model associated with the combined array. The mean squared error for the

estimated process variance for the combined array approach is developed for two cases:

(1) when the model is correctly specified and (2) when the model is underspecified.

Specific combined arrays are evaluated graphically through use of a mean squared

error criterion. Simulation results and additional practical examples outline how

this approach may be used to select an optimal combined array within a particular

experimental situation.











CHAPTER 1
INTRODUCTION



The quality control movement in the 1990s has recognized the need to do more

than find operating conditions that achieve a target value for the mean in order to

improve product quality. Today's industrial statisticians are encountering a need for

developing experimental strategies for achieving some target condition for a charac-

teristic of interest while simultaneously minimizing its variance. To accomplish this,

it is assumed that the experimental factors under consideration separate into a set of

control factors and a set of noise factors. Control factors are those variables under the

direct control of the experimenter. Noise factors are often environmental variables

such as ambient temperature and humidity, properties of raw materials, product ag-

ing, etc. In certain applications they measure how the consumer uses or handles a

product. These noise variables can be controlled by the experimenter during the pilot

plant phase or at the research and development level, but they are difficult to control

at the production level or out in the field.

To motivate this type of optimization problem, let us consider an example dis-

cussed by Box and Jones (1992). Suppose a manufacturer is seeking the best cake

recipe for a cake mix to be sold in a box at the supermarket. Preliminary experi-

mentation has identified three important control factors to be flour (F), shortening

(S), and egg powder (E) and the ranges of F, S, and E that are likely to produce a

good-tasting cake if the mix is baked at a specified oven temperature (T) for a given

length of time (t). Now, it is anticipated that the consumer may have an oven in

which the temperature is biased up or down and that the consumer may overcook







or undercook the cake. Thus, the manufacturer desires a robust recipe so the cake

will taste reasonably good even when the time and temperature of baking, which are

noise variables in this case, differ somewhat from the recommended levels on the box.

Hence, the job of the statistician during the product design stage of the cake mix

is to provide an experimental strategy for finding the settings of the three control

variables that are most robust to changing levels of the noise variables. If we assume

the response will be taste measured on a scale from 1 to 7, with 7 being best, then

we aim to find settings of E, S, and F which achieve a score near 7 with minimum

variability.

The cake mix optimization problem falls under the experimental design method-

ology known as robust parameter design, in which the experimenter attempts to an-

ticipate the maximum amount of noise expected in the process and then finds those

settings of the control variables most robust to the maximum noise. We note that

this robustness to noise variables can only be incorporated into the product design

stage, a fact emphasized by Genichi Taguchi, a pioneer in the methodology for robust

parameter design. Taguchi was among the first practitioners to formalize the idea

that as process variables change, process performance variability changes, and this

change must be incorporated into the statistical design and analysis. In order to im-

prove products and processes, Taguchi recognized that it was no longer sufficient to

concentrate on the mean response alone. In addition, the practitioner must include

product variability as a performance response.

In a first approach to experimentation using both control and noise factors, the

experimental plans introduced by Taguchi have become known as "crossed arrays."

These plans cross one experimental design for the control factors with a separate

design for the noise factors. Crossed arrays have been criticized for requiring a large

number of experimental runs, often wasting resources on the estimation of higher







order control-by-noise (C xN) interactions, and not allowing estimation of control-

by-control (C x C) or noise-by-noise (Nx N) interactions. To address these and other

criticisms and disputes concerning Taguchi's approach, a number of authors have

proposed alternative methods which seek to unite Taguchi's principles with more

traditional statistical design and analysis techniques.

One of the alternative approaches to the construction of a parameter design ex-

periment utilizes a single experimental design in both the control and noise factors.

This approach was introduced by Welch et al. (1990) and discussed further by Shoe-

maker et al. (1991), Lucas (1991), and Box and Jones (1992). Underlying the single

experimental design, now referred to as a "combined array," is a single linear model

in which the response is modeled directly as a function of the control factors and

the noise factors. Shoemaker et al. (1991) noted that the combined array approach

resolves some of the criticisms of Taguchi's crossed array. Specifically, combined ar-

rays can be more economical in terms of run size or number of runs and by proper

construction can allow estimation of all CxC, N xN, and CxN interactions of interest.

In this dissertation we aim to study how the combined array approach can be

used as an experimental strategy for obtaining a target condition on the mean re-

sponse while simultaneously minimizing its variance. The combining of the control

and noise variables into a common array and treating the noise variables as fixed

in the experiment enables us to employ modifications of traditional response surface

methodology. Specifically, optimum control factor settings may be determined by

constructing response surfaces for the process mean and variance. The response sur-

face for the process mean is obtained once we fit a single model in the control and

noise variables. The response surface for the process variance is obtained by applying

the variance operator to the response model in which it is assumed the noise variables

are truly random. Then we substitute the coefficient estimates and denote this esti-

mated response surface for the process variance by f(x). Now the response surface







derived from y(x) is used to locate a contour in the experimental region where the

target is met. Along this contour, we can then obtain the settings in the control

variables that minimize variability. For this purpose we optimize the response surface

derived from (x), the estimated process variance. Since methods for exploring the

response surface generated by the mean of the response and obtaining reliable choices

of control factor settings that optimize (ax) are well documented in the literature, we

will instead focus on T(x), the response surface derived from the estimated process

variance.

This dissertation provides an important contribution to the field of robust parame-

ter design. Combined arrays are currently being utilized by industrial practitioners for

the purpose of designing a quality product or process. This objective is accomplished

by determining settings in the control factors that minimize a measure of loss. Since

this loss is the process variance, a technique is needed for evaluating which combined

array provides "good" estimates of the process variance. However, little consideration

has been devoted to the properties of combined arrays. This research represents the

first attempt at formalizing the properties of the estimated process variance when a

combined array is used. We develop a mean squared error criterion for the estimated

process variance. This criterion will be used to compare the ability of commonly

used designs to support the estimation of r(x) over a region of interest. Thus, if an

experiment is to be run during the product development stage to determine those

settings of the control factors that yield a product that performs well in the presence

of variability, our mean squared error criterion will enable us to recommend to the

practitioner a combined array that minimizes the mean squared error of T(x).

The literature review that follows in Chapter 2 is intended to further acquaint

the reader with the main aspects of Taguchi's methods for robust parameter design

and the criticisms of his methods and to discuss the development of the combined

array approach. In addition, since we will utilize a mean squared error criterion







to compare designs, some background material regarding mean squared error will be

discussed. In Chapter 3 we derive the prediction properties of the combined array and

develop a mean squared error criterion for evaluating combined arrays that are used

to support a model that is presumed to be adequate. In Chapter 4 we do the same for

the case where the fitted model is presumed to be underspecified. In both chapters

we also evaluate combined arrays graphically through the use of a mean squared

error criterion. In Chapter 5 we examine various special cases of the underspecified

model and make recommendations concerning the best combined array to use when

protecting against specific higher order models. Chapter 6 contains conclusions and

potential areas of further research into the combined array and robust parameter

design.












CHAPTER 2
REVIEW OF LITERATURE


2.1 A Review of Taguchi's Parameter Design Methods


The focus of this dissertation is on the properties of the combined array, an ex-

perimental design that was developed to be an alternative to Taguchi's crossed array

design. However, it is important to review some of Taguchi's methods for the de-

sign and analysis of robust parameter design experiments in order to understand the

particular aspects of these methods that prompted the development of the combined

array approach. Kackar (1985) provides a very good review the concept of parameter

design and the quality control methods of Genichi Taguchi. Much of the material

that follows in this section is summarized from Kackar's article.

The evolution of a manufactured product consists of three stages: product design,

process design, and manufacturing. It is the designs of both the product and the

process that play crucial roles in the resulting quality of the product. A good quality

product is one that performs at a target level consistently throughout its life span

and over a wide variety of operating conditions. Define performance variation to be

the amount by that a manufactured product's performance deviates from its target

value under all different operating conditions. Parameter design is an off-line quality

control technique that is used to make a product robust against sources of variation by

identifying settings that minimize performance variation. Parameter design is based

upon classifying the variables that affect a product's performance into two categories:

control factors and noise factors.







The settings of the control factors (design/process parameters) are the product de-

sign characteristics whose nominal values can be easily controlled or manipulated by

the product designer. A vector x' = (x1, ..., xki) of settings of the control variables

defines a product design specification. Noise factors are those that affect product

performance but are difficult or expensive to control. Noise can be classified into two

categories: external and internal. Examples of external sources of noise are ambient

temperature and humidity, suppliers of raw materials, and human variations in oper-

ating the product. Internal sources of noise are those that cause the settings of the

control variables in a manufactured product to deviate from the nominal settings.

The main sources of internal noise are manufacturing imperfections and product de-

terioration. Note that not all sources of noise can be included in a parameter design

experiment because of physical limitations and because of lack of knowledge about

the sources. Let z' = (zi,..., zkn) be the vector of noise variables that are included.

Parameter design experiments are possible because the effects of noise factors may

change with the settings of the control factors. There could be many settings of a

for which the system can perform, on the average, at desired target levels. Among

these, there will be some settings that make the system less sensitive to variation

in the noise variables z. The basic idea in parameter design is to identify, through

exploiting interactions between control variables (C) and noise variables (N), appro-

priate settings of control variables that make the system's performance robust to

uncontrollable variation in z.

A parameter design experiment defined by Taguchi consists of crossing two or-

thogonal arrays: the control (inner) array and the noise (outer) array. Each row

or test run of the matrix for the control variables represents a specific product de-

sign. The rows of the noise matrix represent different combinations of the levels of

the noise factors. Note that by choosing extreme values for the levels of the noise

variables the system is purposely made to be as noisy as possible. The rationale for







this is that settings in the control variables that are determined to perform well in

terms of extreme variability should also perform well under normal conditions. Each

of the n, rows of the control matrix is crossed with the nz rows of the noise matrix

resulting in a crossed array with nxn, rows. For example, suppose there are four con-

trol variables, each having two levels (low = -1, high = 1), and three noise variables

each having two levels (low = -1, high = 1). Attention should be given at this point

regarding the coding of variables. Throughout this dissertation it should be noted

that for convenience all computations will be made assuming that the factors under

study have been transformed from their natural units to design units through the

use of appropriately chosen centering and scaling constants. In this case, we choose

a control array that is 8 x 4 and represents a 24-1 fractional factorial design in the

control variables. We also choose a noise array that is 4 x 3 and represents a 23-1

fractional factorial design in the noise variables. Crossing each of the eight rows in

the control array with the four rows in the noise array yields a 32-run crossed array.

This crossed array is illustrated by Figure 2.1.

In general, Taguchi proposes the use of saturated or nearly saturated two and

three-level experimental designs for both the control and noise arrays. These designs

have been labeled in the literature L4, L8, L9, L12, L16, L18, L20, L24, L27, and

L32 OAs. The 4, 8, 9, etc. refer to the number of design runs, and OA stands for

orthogonal array. Orthogonal arrays are those n x k matrices that have the following

two properties (Rao, 1947):

1. All levels appear an equal number of times in every column.

2. In every pair of columns, all combinations of the levels must occur, and they

must occur an equal number of times.

The notation A(n, k, s, d) has also been used to denote an orthogonal array where n

is the total number of runs, k is the total number of columns, s is the number of levels









X1 X2 3 X4 Z1 Z2 Z3
-1 -1 -1 -1 -1 -1 -1
1 1 -1
1 -1 1
-1 1 1


1 1 -1 -1





1 -1 1 -1





1 -1 -1 1





-1 1 1 -1





-1 1 -1 1





-1 -1 1 1


-1 -1 -1
1 1 -1
1 -1 1
-1 1 1

-1 -1 -1
1 1 -1
1 -1 1
-1 1 1

-1 -1 -1
1 1 -1
1 -1 1
-1 1 1

-1 -1 -1
1 1 -1
1 -1 1
-1 1 1

-1 -1 -1
1 1 -1
1 -1 1
-1 1 1

-1 -1 -1
1 1 -1
1 -1 1
-1 1 1


1 1 1 1 -1 -1 -1
1 1 -1
1 -1 1
-1 1 1


Figure 2.1. Illustration of a Crossed Array







used, and d is the "strength." An orthogonal array of strength d is an n x k matrix

where all sd combinations of levels corresponding to any d factors chosen out of the

k occur an equal number of times. All of Taguchi's designs are orthogonal arrays of

strength 2. A two-level design of strength 2 will allow estimation of the main effects

of the k factors but not pure higher-order terms or interactions between the factors.

Examples of two-level orthogonal arrays of strength 2 are the L4, L8, L12, L16, L20,

L24, and L32 designs. The L4, L8, L16, and L32 designs are the familiar 22, 23, 24,

and 2s resolution III fractions of 2k factorial designs. The L12, L20, and L24 designs

are two-level Plackett-Burman designs. The L9, L18, and L27 designs are three-level

Plackett-Burman designs. Note that these three-level orthogonal arrays of strength 2

will allow estimation of main effects and pure quadratic terms but not interactions.

Lastly, note that Taguchi's experimental designs are formed by crossing the control

array with the noise array, and this crossed design is such that all first order and

higher C x N interactions are estimable. The crossed array displayed in Figure 2.1

is an L8 control array crossed with an L4 noise array.

The analysis of a crossed array uses the nz values of the response variable (or qual-

ity characteristic) obtained for each of the n, rows of the control matrix to compute

n, values of a performance statistic that Taguchi calls a signal-to-noise (SN) ratio.

A performance statistic is a function of the control variables that estimates the effect

of noise factors on the characteristic of interest. The values of this statistic are then

used to predict optimal settings of the control variables.

The identification of optimal control factor settings requires specification of a

criterion to be optimized. One such criterion is the expected loss. Let Y be a

value of the performance characteristic of interest and suppose its target value is

t. Further, assume that Y is a function of the control and noise variables, i.e.,

Y = f(x, z). The goal of parameter design is to choose the setting of x that minimizes

the expected loss caused by deviations of Y from the target value t. This expected







loss is given by R(x) = EzL(Y, t) where L is a loss function. Taguchi recommends

that the loss be measured by the quadratic loss function L(Y, t) = K(Y t)2 where

K is some constant. In practice, however, Taguchi maximizes SN ratios instead of

minimizing R(x) and generally gives no connection between these two optimization

problems. Leon et al. (1987) show that for certain underlying models for the response,

maximization of the SN ratio does indeed correspond to minimization of the expected

quadratic loss.

The forms of the SN ratios suggested by Taguchi depend on whether "target-is-

best," "smaller-is-better," or "larger-is-better." The SN ratio for the target-is-best

case takes into account the relationship between the mean and the variance. If the

variance increases linearly with the mean, then the appropriate SN ratio is given by

SN, = 101oglo(y72/s)

where y, and s5 are the sample average and sample variance of the nz response values

at the ith control setting, i = 1,2,..., n If the mean and variance are independent,

then the appropriate SN ratio is

SN, = logo s,.

For the cases of smaller-is-better and larger-is-better, the SN ratios are respectively

given by

n, n,z /y
SN= -101oglo(Ey /nz) and SN,=-101og0o--('-).
j=1 j=1 nz

In each case described above it is desired to maximize the appropriate SN ratio.

To carry out this optimization, Taguchi introduces the concept of adjustment factors.

These are special control variables whose values can easily be changed to reduce an

adjustable component of the total variation about the target without affecting the

nonadjustable component of the total variation. In other words the adjustment factors







are those control variables that only affect the mean. The remaining control factors

affect both the mean and variance or just the variance. The division of control factors

into two groups x = (d, a) enables the analyst to minimize the variability using the

main control factors d and then adjust the mean to its appropriate target using

the adjustment factors a. In summary, Taguchi recommends the following two-step

optimization procedure:

Step 1: Find the setting d = d* that maximizes the SN ratio.

Step 2: Adjust a to a* while d is fixed at d* in order to minimize the bias.

Leon et al. (1987) provide the link between maximizing the SN ratio and minimiz-

ing the loss R(x). They show that if Y = f(x, z) has the multiplicative form given

by Y = y(d,a)e(z,d) where E[Y] = p(d,a) is a strictly monotone function of each

component of a for each d, then the use of Taguchi's SN ratio in a two-step procedure

leads to minimization of quadratic loss. To see this, the authors note that the mul-

tiplicative model implies that the squared coefficient of variation, var[Y]/E2[Y], does

not depend on a (or approximately, for log Y, a affects location but not dispersion).

Also, to link the SN ratio with quadratic loss R(x), the authors argue as follows:

1. First note that we can find (d*, a*) such that R(d*, a*) = minda R(d, a).

2. Define the general two-step optimization procedure:

Step 1: Find d* that minimizes P(d) = mina R(d, a).

Step 2: Find a* that minimizes R(d*, a*).

3. If we use quadratic loss and the multiplicative model, then P(d) is a decreas-

ing function of Taguchi's SN ratio and thus the above two-step procedure is

equivalent to Taguchi's.

Taguchi's method for conducting a robust parameter design experiment and sub-

sequently maximizing a signal-to-noise ratio can be summarized as follows:







1. Identify initial and competing settings of the control variables, dividing them

into adjustment and nonadjustment variables if possible. Identify the noise

variables and specify the ranges of the noise variables within which the product

performance is desired to be insensitive.

2. Use orthogonal arrays to construct the control and noise variable matrices.

3. Conduct the parameter design experiment and evaluate the performance statis-

tic (SN ratio) for each test run of the control matrix.

4. Use the computed values of the SN ratio to predict new settings of the control

variables that maximize the SN ratio.

5. Conduct a confirmatory run to check that the new settings do indeed improve

the performance statistic.

Lastly, note that the foregoing discussion was concerned with the product design

stage in the evolution of a manufactured product. In the process design context,

parameter design experiments are used to identify settings of controllable process

variables at which the effect of manufacturing variation is minimum. The control

variables are thus process variables whose operating standard can be chosen by the

process designer, and manufacturing variations are sources of noise.

2.2 The Combined Array: An Alternative to Taguchi's Crossed Array


Recall that Taguchi's approach uses two experimental designs: one for the con-

trol factors and one for the noise factors; and for each row in the control array the

"replicate" observations generated by the noise array are reduced to a loss statistic

or related SN ratio that can be used to compare specific product designs defined by

the settings of the control factors. The observed losses are then modeled. A common

criticism of this approach is that the model often includes only the main effects of







the control variables and little attention is paid to the effect of interactions among

the control variables on the underlying quality characteristic Y.

Welch et al. (1990) proposed a new approach to parameter design that models

a quality characteristic Y(x, z) as a function, typically polynomial, of both x and

z. The single experimental design for both types of parameters generally requires

far fewer runs than crossed arrays. The method of Welch et al. (1990) involves the

following six steps:

1. Design a single experiment to predict the response of interest as a function of

the control factors x and the noise factors z, where the response model is given

by





The fjs are assumed known functions of x and z and the /3js are unknown

model coefficients to be estimated from the data. Hence the combined array

approach selects designs on the basis of their ability to estimate the appropri-

ate coefficients in the presumed model for the response. Attention should be

given at this point regarding the application in this article. The authors used

computer simulation to model a clock skew circuit design and thus the error

term above represents systematic departure from the assumed model because

the computer response Y is deterministic.

2. Predict Y(x, z) using the ordinary least squares estimator


Y(x,z) = E f(x, z).


3. Once the terms in the model are appropriately estimated, the experimenter can

predict the value of the loss statistic R(x) for a given x. For the example in

which a clock skew circuit design is modeled using computer simulation, the







target for the simulation output Y(x, z) is zero, and following Taguchi, the

expected squared error loss is given by


R() = JY2(mx,z)g(z)dz,


which averages a measure of loss over the distribution g(z) of the noise factors.

Denote the predicted value of this loss statistic by


R(x) fY2(x, z)g(z)dz.


Welch et al. (1990), for the purposes of their example, assume a uniform distri-

bution for the noise factors.

4. Minimize the estimated expected loss !R(x) as a function of x.

5. Conduct a confirmatory experiment to evaluate R(x) by fixing x at the values)

found in step four and varying z according to g(z).

6. Iterate if necessary. The optimization step may suggest for example that the

design region needs to be shifted.


Shoemaker et al. (1991) suggest complementing the formal analysis of the expected

loss in step four above with graphical data-analytic techniques. Specifically, the

authors recommend the use of control-by-noise interaction plots to help determine

the control factor settings that dampen the effects of the individual noise factors.

To illustrate and compare their method with others, Welch et al. (1990) pro-

vided an example in which the goal was to determine the widths wl,..., w6 for the

transistors that gave the smallest clock skews (si and s2) in the presence of process

variability. Three strategies were employed. The first was a 60-run experiment for the

skew models from which two specified loss statistics were predicted indirectly (this is







the new method of Welch et al., 1990). The second was a (40 x 5)-run crossed array

(i.e., separate control and noise array as advocated by Taguchi) for modeling the loss

statistics directly. The third was a hybrid strategy that used the (40 x 5) experiment

from the second strategy but modeled the skews to predict the losses indirectly as in

the first strategy.

A comparison of the predicted values of the loss statistics with the actual values

obtained in confirmatory runs for the best three vectors of control factor settings

resulted in the following conclusions:

1. Strategy one predicted the skews accurately enough to give predictions very

close to the actual losses determined by confirmatory experiments. The pre-

dictions were more reliable and the actual losses were smaller than in strategy

two.

2. Strategy one gave virtually identical results as strategy three, indicating that

small experiments can be adequate.

3. Modeling the skews (i.e., the response) appeared to be more accurate than

modeling the loss statistics directly.

4. Strategy one gave the same best-three circuit designs for the two loss statistics

considered.

5. The authors point out that choosing models for complex loss functions, as in the

second strategy, is difficult and often not very intuitive. When the response was

modeled directly, certain control factor interactions were justifiably omitted, but

the same interactions may not be negligible when modeling the loss statistics.

A full second-order model in the widths was chosen for both loss statistics.

Also, modeling the response generally admits simpler models. Furthermore,

collapsing the data to a loss statistic could hide important relationships in the

data.







The method of Welch et al. (1990) can be summarized simply as follows: combine

control and noise factors into a single array; model the response rather than expected

loss; and approximate a prediction model for the loss based on the fitted response

model. The single array for the control and noise factors has since become known as

the combined array.

2.3 Run Savings and Estimation Flexibility of the Combined Array


Shoemaker et al. (1991) further developed the response-model/combined-array

approach suggested by Welch et al. (1990). The authors' contributions to parameter

design are motivated by what they consider to be the following disadvantages of

Taguchi's loss model approach and corresponding crossed array designs:


1. A very large number of runs may be required because the noise array (NA) is

repeated for every row in the control array (CA).

2. Many degrees of freedom (df) are dedicated to estimation of all interactions

between control factors and noise factors.

3. The focus is on modeling the loss, or related SN ratio, which is often a nonlinear,

many-to-one transformation of the response (i.e., performance characteristic) Y.

Even when Y follows a linear model in the control and noise factors, it is less

likely that the loss can be modeled well by a low order linear model even if data

transformations are employed.

4. Modeling the loss may also hide some of the relationships between individual

control and noise factors.


To elaborate on the second disadvantage stated above, Shoemaker et al. (1991)

show that for the general crossed array case, the effects that are estimable in a CA x

NA design correspond to those that are estimable in CA, those that are estimable







in NA, and the generalized interactions between the estimable effects in the control

array and the estimable effects in the noise array. This fact illustrates one of the

limitations of the crossed array approach. If the control array has n., runs and the

noise array has n. runs, then the CA x NA design has nnz runs, and of the nnz 1

df in CA x NA, (n, 1)(nz 1) df are used for estimating two-factor and higher

order control-by-noise (C x N) interactions. In practice, however, some of the C x N

interactions may be judged a priori to be unimportant; so, the large run size of the

crossed array is often unnecessary. Also, there is no flexibility to reallocate some of

the df to estimate control-by-control (C x C) interactions that may be important.

Taguchi's three-level designs are also unable to estimate (C x C) interactions.

The inflexibility and large run size of the crossed array can be avoided by combin-

ing the control and noise factors in a single array. This allows for a smaller experiment

since a more highly fractionated design may be chosen. Also, with suitable choice

of a defining relation, a fractional factorial design can be constructed in which some

C x C interactions will be estimable if certain C x N interactions are assumed neg-

ligible. To illustrate this, Shoemaker et al. (1991) provide an example to compare

the estimation flexibility of a crossed array design with a combined array design.

Three control factors and three noise factors, each with two levels, were used. To

construct the crossed array, a 2R'1 CA was crossed with a 231' NA. The resulting

16 run, six variable crossed array is equivalent to a 26-2 fractional factorial design.

On the other hand, by combining the six factors into a single array, a 2'62 fractional

factorial design can be constructed using defining relations that allow estimation of

certain C x C interactions. From the example, it was observed that the combined

array is superior to the crossed array because it allows more effects to be estimated

under weaker assumptions.

In conclusion, Shoemaker et al. (1991) illustrate that the run savings allowed by

the combined array format comes from the flexibility to estimate those effects that







are most likely to be important, rather than forced estimation of a large number of

control-by-noise interactions.

2.4 Analysis of a Combined Array Experiment


The development of optimization strategies for robust parameter design appli-

cations is important for the statistician in industrial practice. Recall that robust

parameter design experiments are possible because the effects of noise factors may

change with the settings of the control factors. Since the variance fluctuates over

the region of interest, Taguchi was the first to emphasize the possibility of modeling

the variance. The basic idea in parameter design is to identify, through exploiting

interactions between control variables (C) and noise variables (N), appropriate set-

tings of control variables that make the product or process performance robust to

uncontrollable variation in z. In other words, although there may be many settings

of x for which the product or process can perform, on the average, at desired target

levels, the goal is to determine the setting that makes the product or process less

sensitive to variation in the noise variables z, and this is done by exploiting C x N

interactions.

Myers, Khuri, and Vining (1992) show that if the model for the response contains

C x N interactions, then the process variance is a function of the levels of the con-

trol variables; on the other hand, if C x N interactions are not included (i.e., they

are deemed to be unnecessary or insignificant), then the process variance does not

depend on x, and in this case the "robust" parameter design problem is nonexistent.

Consequently, when interactions between control and noise variables are found, the

process variability depends on the specific values of the control variables through the

C x N interactions included in the underlying model. The interaction structure then

can provide a reasonable basis for estimating the process variance where the noise

factors are truly random variables.







As an example, suppose that for a particular combined array experiment, there

is one control variable and one noise variable. If there is no interaction between the

control and noise variables, then the fixed effects model for the experiment is given

by

y= -o + xi 71z +yii c-.

Now, in reality, z1 is a random variable with a variance ao, and thus the process

variance is given by

var(y) 1= 72( + C2.

Notice here that the variance does not depend on the design levels in x\. Hence there

is no need to incorporate the process variance as a part of the performance criterion

for parameter design. On the other hand, consider the fixed effects model for the

combined array experiment where it is assumed that there is interaction between the

control variable and noise variable. This model is given by


y = /3o + iXi + 71zi + Aixxzi + e. (2.1)

Again assume that in reality the levels of the noise variable are random in the process

so the process variance is given by


var(y) = (71 + Aixi)2 a + a (2.2)


The process variance is now a function of the levels of the control variable, and thus it

is important that the process variance be involved in the response surface performance

criterion.

The example above has provided motivation for including the process variance

in the choice of optimum settings for the control variables. Note that the variance

function (2.2) is essentially a quadratic response surface. The underlying model (2.1),

where the noise variables are treated as fixed effects, provides a basis for estimating

'71 and A,, and thus allows estimation of the process variance. The estimated variance







response surface is given by

var(y)= (71+ +iXx)2 + ^2

where we assume that or is known. Now, assume that E(e) = 0, and without loss

of generality assume that E(zi) = 0. Then E(y) = 3o + A/1x, and so the estimated

response surface for the mean is given by

(x) = 3o + /xi.

As a result, optimum conditions may be obtained through simultaneous exploration

of these two response surfaces. One method for exploring the response surfaces is

through the construction of contour plots. Joint consideration of these contour plots

can help determine the optimum conditions where the target value for the response

is obtained with minimum variability. More formal methods for dual response opti-

mization can be found in Myers and Carter (1973), Vining and Myers (1990), and

Del Castillo and Montgomery (1993).

2.5 The Mean Squared Error Criterion


Let us suppose that the true functional relationship between a response r7 and k

continuous variables x1, x2,... Xk is given by yf = g(x, 0). When the form of the true

functional relationship is unknown, the objective is to approximate within a given

region R of the k-dimensional space of the variables, the function g(x, 0) by some

function f(a, /3). Since the fitted model may be an inadequate representation (i.e., an

underestimate) of the true response surface, there is a need to consider two sources

of error. Error can arise both from sampling variance and bias introduced by model

inadequacy. In other words, we explain the discrepancy between the true function

r(x) and the fitted function y(x) as being due to "variance error" and "bias error."

In this problem setting, Box and Draper (1959, 1963) propose that the design D,

i.e., an n x k matrix whose (ui)th element, xu,, is the value of the ith input variable







used in the uth experimental run (u = 1, 2,..., n; i = 1,2,... k), should be chosen so
as to minimize the quantity

nK
J- JR E[y(x) -r(x)]2dx


where dx = dxedx2 ... dxk, K = [fR dx]-', E[ (x) i(x)]2 is the mean squared error
of the predictor y(x), o2 is the experimental error variance, and n is the number of
design runs. Here we interpret J as being the mean squared deviation of y(x) from
the true response 7r(x), averaged over the region R and normalized with respect to
the number of designs runs and the experimental error variance.
It can easily be shown that J separates into contributions from the "variance
error" and the "bias error" according to


J=nK- var[ (x)]dx + ))2d= V + B.


This sum of the average variance and the average squared bias over the region R has
become known as the Box-Draper "J-Criterion" or the integrated mean squared error
(IMSE) criterion. The utility of the J-Criterion is that it represents a compromise
between purely variance based or purely bias based criteria that had dominated the
literature on the construction of optimal designs prior to Box and Draper's work.
Box and Draper (1959, 1963) apply the J-Criterion to the setting in which f(x, 3)
is assumed to be a polynomial of degree d, while the true functional form g(x, 0) is
a polynomial of degree d2 = d1 + 1 where di was taken to be equal to 1 and 2.
Designs having minimum integrated mean squared error were obtained for various
ratios of V/B and compared with "all variance" (B = 0) and "all bias" (V = 0)
designs. The designs were constructed to satisfy certain moment conditions derived
from the minimization of J. It was concluded that if the practitioner is concerned
with protecting against model inadequacy, then the design that minimizes J is close
to the "all bias" design.







In the combined array approach, recall that we assume a model for the response

and then approximate a prediction model for the loss based on the fitted response

model. In practical circumstances, it is reasonable to believe that the assumed model

is always to some extent incorrect. We will develop a J-Criterion for the estimated

process variance, in light of Box and Draper's work, not for constructing a new class

of optimal designs, but rather as a criterion for choosing between already existing

designs that have been used for constructing combined arrays. Again we emphasize

that in the robust parameter design problem that we are studying, the focus is not on

minimizing the bias and variance of the predicted response y(x). Rather, the focus

is on the estimated process variance, var(y), which we denote by (ax). The rationale

behind this is that we assume that a set of operating conditions for which the mean

response achieves a target value has been found, and the objective is to determine

from this set the setting x that achieves the target value with minimum variability.

These conditions are such that (a) is minimized. Hence, we focus on the properties

of N(a).

It is hoped that this review has set the stage for the following chapters in which the

reader will be shown how the abundance of literature in traditional response surface

methodology can be utilized and extended to develop properties of the combined

array. Once these properties are established, the theory will be laid for the industrial

practitioner who must choose the most appropriate experimental design to run for

his or her particular robust parameter design problem.












CHAPTER 3
ANALYSIS OF THE COMBINED ARRAY
WHEN THE ASSUMED MODEL IS CORRECT


3.1 Introduction


In the combined array approach designs are selected on the basis of their abil-

ity to support the estimation of the coefficients in the presumed model for the re-

sponse. Once these terms are estimated, the analyst is able to estimate the appropri-

ate squared error loss. Typically, an additive error structure is assumed, in which case

the process variance is the appropriate loss function. In a robust parameter design

problem the predicted response and the predicted variability are equally important

for determining the optimum control factor settings. We assume the experiment will

be run in a region where the predicted response is likely to achieve a target value, and

our goal is to determine where in the region the variability around the target value is

minimized. Hence the estimated process variance becomes the performance criterion

that we focus on.

In this chapter, we will determine the properties of the estimated process variance

under the assumption that the fitted model is correct. The statistician would like

to choose an experimental design that will produce accurate and reliable estimates

of the process variance over the entire experimental region. Hence, we will develop

a mean squared error criterion in order to evaluate the properties of the estimated

process variance obtained from a given design. For a given model, a series of sim-

ulations will be conducted to illustrate the use of the mean squared error criterion

for comparing commonly used designs to fit the given model. An additional practical







example outlines how this approach can be used to select an optimal combined array

among a given set of experimental designs. The mean squared errors will be displayed

graphically in order to make conclusions and recommendations regarding competitor

designs.

3.2 Properties of the Estimated Process Variance


For a combined array experiment suppose that the model for the process is given

by

y = 3o + x', + z'y + x'Az + e (3.1)

where

/o is the intercept,

x is a p x 1 vector representing the appropriate polynomial expansion of the

control variables for the assumed model,

3 is the p x 1 vector of unknown parameters associated with the control factors,

z is a q x 1 vector representing the appropriate polynomial expansion of the

noise variables for the assumed model,

7 is the q x 1 vector of unknown parameters associated with the noise factors,

A is the p x q matrix of unknown parameters associated with the control-by-

noise interactions (some elements of A may be zero), and

e is the random error.

We shall assume that the es are independent and identically distributed with mean

zero and variance ucy. For the purposes of the particular planned experiment, we

shall treat the noise variables as having fixed effects, and thus there is no distinction

between the roles of the control and noise variables. In the process, however, the







noise factors are truly random variables with E(z) = p and var(z) = V, where we

shall assume V is known. Furthermore, the error terms and noise variables shall be

assumed to be independent.

Assume the appropriate combined array experiment has been run to estimate the

linear model given by equation (3.1). The analysis focuses on both the response and an

appropriate loss function. Leon et al. (1987) show that if the goal of experimentation

is to find conditions for which a target value is achieved with minimum variability

and if an additive model as in (3.1) is assumed, then the appropriate loss function

is the process variance, which we shall denote by r(x). Since the noise variables are

considered to be random in the process, the value of the process variance at a point

x inside a region of interest R is, for model (3.1), given by

T(x) = var(y) = var[(Y' + x'A)z + e] = [7' + x'A]V [7' + x'A]' + a2

where V is the q x q variance-covariance matrix for z. Now if we define O(x) =

7 + A'a, we can write

r(x) = var(y) = 0'(x)V (ax) + a~l. (3.2)

Since O(x) is unknown, we use the results of the combined array experiment to

estimate O(x) with 4(x). The estimated process variance is then given by f(x) =

' (x)V<>(xs) + <2.. The aim is to find control factor settings x that minimize (ax),

subject to the constraint that the predicted response is at the target value.

The derivation of prediction properties of f(x) will make use of the following

more convenient notation. Let A' = [y, A'] and x'a = [1,x']. Then O(x) may be

reexpressed as O(x) = A'xa. Also for convenience, let A = vec(Aa), where the

vector of the matrix Aa is the q(p + 1) x 1 column vector formed by stacking the

columns of Aa one under the other. Using this notation we can rewrite the model

given by (3.1) in the more general form









y = P +E

where

y is the n x 1 vector of responses,

P is the n x m model matrix (m = (p + 1)(q + 1)) whose elements are known
functions of the settings of k = kc + k, input variables determined by an n x k
design matrix, i.e., P = [1, X, Z, XZ] where X and Z are n x p and n x q
design matrices in the control and noise variables, respectively, and XZ is an
n x pq matrix of cross products,

P' = [/o, /3', A'] is the m x 1 vector of unknown parameters, and

e is the n x 1 random error vector.

If we assume that E[E] = 0 and var[e] = olI and use ordinary least squares (OLS)
to estimate the model coefficients, the resulting estimate of the coefficient vector is

, = (P'P)-1'Py. By the properties of OLS estimation and under the assumptions
that the model is correct and the noise variables have fixed effects in the experiment,

E[ ] = 4 and var[ ] = cy(P'P)-'. Note that /(x) = ,+A 'a is a linear combination

of the elements of ,, and so E[ (x)] = O(). Thus,

E[(x)] = E[k'(x)V(x) + ]

= tr(VE) + 0'(x)V((x) + cU (3.3)

where E = var[4(x)]. From (3.2) and (3.3) we note that

bias[f(xa)] = E[f(ax)] r(x) = tr(VE). (3.4)

Note that V and E are both variance-covariance matrices. Consequently, V and
E are positive definite. Thus, it follows that tr(VE) > 0 and so bias[(xa)] > 0







(Graybill, 1983, p. 397). We thus obtain an unexpected result, which is that ?(ax) is

a biased estimator of r(x) even though the presumed model is believed to be correct.

Next if we assume the random error vector e is distributed as multivariate normal

with mean 0 and variance a2I then the variance of f(xa) is given by



22
var[(a)] = var[ck (x)Vo(x) + &,]

= 2tr(VE)2 + 40'(x)VYEV(ax) + 2-u (3.5)
n- m

where n is the total number of design runs and m is the number of unknown coefficients

in the model.

3.3 Expressions for the Bias and Variance of the Estimated Process Variance


To derive a more explicit expression for the bias, assume the correct model is

given by (3.1). From least squares estimation of the equivalent form of the model

given by y = Pi + e, we noted above that ^ = [#o, 3', A']' = (P'P)-lP'y with

var[^] = a2(P'P)-'. Let C be the appropriate q(p + 1) x q(p + 1) submatrix of

a7(P'P)-1 corresponding to the variance-covariance structure of A. Note that C is

symmetric and can be partitioned into q2 blocks of (p + 1) x (p + 1) matrices of the

form
C11 C12 CIq' q
C 2 C12 C22 C2q (3.6)

C.q 2q C qq

Since (ac) = A'am, and A = vec(Aa) it follows that (ax) = (Iq 0 X'a)A, where

denotes the Kronecker cross product. Hence


E = var[(B(x)] = (Iq 0 X')C(Iq 0 )' = (Iq 0 X'a)C(Iq 0 a),







where we used the general result that (A B)' = A' 0 B'. It is easy to see that

S {aYxaCijXa} for i,j = 1,2,...,q and let V = {ciij} for i,j,= 1,2,...,q. Now

recall that for two q x q matrices A and B, tr(AB) = Zf Z. aijbji. Thus, it

follows that

q q
bias[(.(x)] = tr(VE) = a 2 orijsXaCjiXa
i=1 j=l

q q-1 q
T 7, u ( Cjx.+ 2 2 ^ C3 j a (3.7)
i=1 i

where we used the fact that aij = Oji and Cj = C'.. Further simplification of the

bias expression depends upon the form of the Cij matrices which in turn depends on

the model matrix P for a given design.

Next, consider the variance of f(x). Utilizing the result that O(x) = (Iq 0 x')A,

equation (3.5) can be expressed as

274
var[( (x)] = 2tr(VE)2 + 4A'[Iq 0 Xa]VEV[Iq 0 x'a]A + 2- (3.8)


Our goal is to evaluate var[f(x)] for any given xa, but observe that var[-(x)] is a

function not only of the particular location in the design region but the unknown

parameter vector A as well, which presents obvious difficulties. Our approaches for

overcoming these difficulties will make extensive use of the Euclidean norm of the

vector A, which we shall denote by r = JAI. We shall assume that we have some

reasonable idea about the value of K. Next we use a result proven by Stroud (1971).

This result states that if U is the surface of the unit sphere and ds8 a differential on

this surface, then the average value of the quadratic form A'AA over all possible A's

with length K is given by

K2 fu A'AAdsA (A
dtr(A)
fu ds \ nN







where n\ is the number of nonzero elements of A. Applying this result to equation

(3.8), the average value for var[,a(x)], denoted by varA[f (X)], over all possible As with

length n, is


)] fu var[ (x)]dsA
f dvarA


= 2tr(VVE)2 + 4 tr[(Iq 0 a)V V(I, 0 X')] + 2u4 (3.9)
nx n m

where nA is the number of elements in A not presumed a prior to be 0. From (3.4)

and (3.9) we can obtain

MSE[f(x)] = [bias(f(x))]2 + varA[(X)]

4c2
= [tr(VE)]2 + 2tr(VE)2 + --tr[(Iq xa)VEV(Iq 0 x')]
nA

2,4
+ -- (3.10)


which can be readily calculated at any point Xa.

3.4 Comparing Designs Based on the Mean Squared Error Criterion


A reasonable method for comparing designs is to evaluate the average or integrated
mean squared error (IMSE) of f(x) for each design. The IMSE is given by


IMSE[f(a)] = K I MSE[f(x)]dx


= K J[tr(VE)]2ddx

f 4r2 2U4
+ K [2tr(V )2 + tr[(q )Vx V(Iq, 0 ')] + ]d (3.11)
R nA n m

where R is the appropriate region of interest and K = [fR dx]-1. Thus, analogous to

the Box-Draper (1959) J = B + V criterion, we note that the average mean squared







error of (xa) is the sum of the average squared bias of (ax) and the average variance
of (x), where average means averaged over the region R. The IMSE can be found

using Monte Carlo integration techniques.
A series of simulations will be conducted to illustrate the proposed IMSE criterion
for comparing designs. We will consider a robust parameter design experiment con-

sisting of three control variables and two noise variables. Suppose we use the model

given by (3.1), and assume the model is first order in the control and noise variables

and contains first order C x N interactions. This model, in expanded form, is given by

y = /o + /lXl + 4/2x2 -+ /3X3 + 71'YZl + 7'22 + Alx1lz1 + A21X2Z1 + A31X3Z1

+ A12xiZ2 + A22X2Z2 + A32X3z2 + e. (3.12)

From equation (3.11), it can be noted that computation of the integrated mean

squared error of a(x) at any point x depends on V, the assumed known variance-

covariance matrix for the noise variables, and the value of K2. For this example with

two noise variables, let orn, o12, and 022 denote the elements of the V. Since our

method assumes the noise variables follow a normal distribution, it seems reason-

able to estimate all and a22 using (Rfi/4)2, where Ri denotes the estimated range

(specified by the experimenter) of the it noise variable, i = 1, 2. If the noise vari-

ables are each studied at two levels in the experiment, conveniently denoted by -1

and +1 in coded (i.e., centered and scaled) design units, then 1, = 2, i = 1, 2, and

aO' = 622 = (2/4)2 = 0.25. (Recall that throughout this dissertation all computa-

tions are made assuming that the factors under study have been transformed from

their natural units to design units through the use of appropriately chosen center-

ing and scaling constants. We prove in Appendix A that the mean squared error

of (x) is invariant to our choice of scaling in the noise variables.) Next, since

012 = cov(zl,z2) = p12(0o11622)2 = P12/4, where p12 is Pearson's coefficient of corre-
lation, we have transformed the problem so that the experimenter now only needs to







have knowledge of the correlation structure of the noise variables. In other words,
the matrix V can be expressed as a constant times the correlation matrix:


V = 1 1 p12
4 p12 1 J


To investigate the effect of the correlation structure of the noise factors on the bias and
mean squared error, we shall consider a sequence of values P12 between -1 and 1. Next,

note that if iR? is misspecified, then it can easily be shown that var[zi] = l(Ri/Ri)2.
Hence, to examine the impact of misspecifying the range on the bias and mean squared
error, we shall consider three cases:

1.
2. II = 1; 722 = .25 < (Ri = R1/2,R2 = R2);

3. oai = 2.25; 022 = .25 = (Ri = R1/3, R2 = R2)

In cases two and three the corresponding V matrices are respectively given by


V = 1 [2 P12 and V= [3 3 P12
2 p12 1/2 and 4 p12 1/3

Note that we need not consider the possibility of misspecifying the range of both

noise variables by the same constant k. Clearly, if R1 = kR1 and R2 = kR2 then
V = var[z] will be multiplied by k and the resulting mean squared error in (3.10)

will be multiplied by k2. Thus, although the magnitude of the mean squared error
will change, the relative ranking of the competitor designs with respect to the mean
squared error will not change.
Corresponding to each of the above three cases, we shall consider two values of
K2 = 1Al2. In this example A' = (71, A, A21, A31i, 72, A12, A22, A32) where the elements

of A are the coefficients of the coded noise factors and control-by-noise interactions.







Without loss of generality assume that a, = 1. Now if we assume that the magnitude

of each element of A is four times the magnitude of ao, then n2 = na(16au) =

(8)(16oa) = 128. Secondly, if we assume that the magnitude of each element of

A is twelve times the magnitude of ao, then c2 = (8)(144U2) = 1152. Consequently,

the constant 402/nx in equation (3.11) assumes the values 64 and 576. Thus, the three

V matrices and two values of K2 represent a total of six different cases of interest.

These six cases intentionally cover a wide variety of states of nature that will allow

for a fair comparison of different combined arrays. In addition, this will enable us to

investigate whether the methodology developed in this paper for comparing designs

is "robust" to the various conditions that the practitioner may come across.

To illustrate the use of the IMSE criterion, we will compare four combined arrays

that support the estimation of model (3.12):

Design 1: 16-run Resolution V Fractional Factorial Design

Design 2: 20-run Plackett-Burman Design

Design 3: 24-run Plackett-Burman Design

Design 4: 24-run Orthogonal Main-Effect Plan.

These designs are easy to construct, and for strict first-order models they are or-

thogonal. Since model (3.12) contains interactions between the control and noise

variables, designs 2-4 are nonorthogonal, but they still support the estimation of all

of the parameters in the model.

3.4.1 Design 1: 16-run Resolution V Fractional Factorial Design


A full 2k factorial design, expressed in design variables, is every possible combina-

tion of 1's for the k factors of interest. The full 2k factorial design allows estimation







of all main effects and first order and higher interactions among the variables. How-

ever, when fitting a model such as (3.12), we are not interested in estimating higher

order interactions; so a fraction of a 25 factorial design may be used. A Resolution

V fraction of a factorial design is an orthogonal array where main effects are aliased

with four-factor interactions, and first order interactions are aliased with three-factor

interactions. By "aliased" we mean that the main effects cannot be estimated inde-

pendently of the four-factor interactions, and the first order interactions cannot be

estimated independently of the three-factor interactions. However, if model (3.12) is

correct, these higher order interactions are assumed to be negligible. For the five-

variable example that we are considering, the 25-1 Resolution V design matrix is

illustrated in Figure 3.1.


x1 X2 X3 zl z2
1 -1 -1 -1 -1
-1 1 -1 -1 -1
-1 -1 1 -1 -1
-1 -1 -1 1 -1
-1 -1 -1 -1 1
1 1 1 -1 -1
1 1 -1 1 -1
1 1 -1 -1 1
1 -1 1 1 -1
1 -1 1 -1 1
1 -1 -1 1 1
-1 1 1 1 -1
-1 1 1 -1 1
-1 1 -1 1 1
-1 -1 1 1 1
1 1 1 1 1

Figure 3.1. 16-run Resolution V Fractional


Factorial Design


Computation of the mean squared error of r(ax) given by equation (3.10) is sim-

plified for the case when the combined array is a Resolution V or higher, two-level

fractional factorial design. For such a fractional factorial design that is appropriately







constructed to estimate the C x N interactions in the model, the matrix C given by
equation (3.6) has elements Cij = 0 for i 7 j and Ci = 4I. Consequently, the bias

of T(s) given by equation (3.7) reduces to

2 q
bias [f(a)] = --' E cr
n i=1

where r2 = (x'xa,) = 1 + EZ1 xz. Thus, note that the bias depends the variance of
the zi's, but not on the correlation between zi and zj. In other words, the bias is

constant over values of pij. Under case 1 in the example above, where p = kc = 3,

q = kn = 2, orn = var(zi) = .25, 0C22 = var(2) = .25, n = 16, and o2 = 1, we have


bias[f(x)] = 6(.25 + .25)= .


Since Xa = (1, x1,x2, X3)', it follows that r2 = 1 + X2 + x2 + x2. Thus, the average

bias can be determined as follows:


K bias[V(a)]dx = K2 [(1 + x + x2 + x)d

S1l f1l f-11(1 + + x2 + x)dxidx2dx3
32 f l- fi fdxldx2dx3

24 1
2 ) = = .0625.
32(23) 16

Similarly, it can be shown that the average squared bias portion of the IMSE is equal

to 0.0041667. Next consider the maximum value of bias[a(x)]. In general, since we

assume the control factors have been coded such that -1 xi : 1 i = 1,..., k, the

maximum value of r2 is kc + 1, and so the maximum bias is


S(kc + 1) q
n i=1







For the example being discussed it is easy to show that the maximum bias is equal

to 0.125. It is not surprising to note that the bias is maximized at the perimeter of

the region of interest. Also note that in general, the bias of a(x) can be large if the

design is highly fractionated, if the inherent process variability o-2 is large, or if there

is substantial variability in the noise variables.

Next consider the variance of "(x). For the Resolution V fractional factorial

design, recall that the matrix C given by equation (3.6) simplifies to -2I, and thus


E = var[4(x)] = (Iq 0 X')C(Iq 0 Xa)

2
n
2 2 2
X1 X I
n n

Thus varA[T(x)] given by equation (3.9) reduces to


294r4 4 + 2 2 2 4
vaA[r(x)] = r tr(V2) + 4 -tr[(I, x)V2(I, )] + ---
n nA na n-m


24r 442 2 2 224
= 2ortr(V2) + 4r, tr [V2 0 XaX'a +
n nA n n- m

S2o4r4 V2) + 4K2 o2r2 (V2,)( + 2o74


= ( + )r tr(V2) + ( n-
n- nn + (- m)

2 4 2-4
n n, n i=1j=1 (n -m)


where we used the general results that tr(xa'aX) = r2, tr(A 0 B) = tr(A)tr(B),

and tr(V2) a= =1 ,=J oj for a q x q symmetric matrix V having elements aij,

i,j = 1,2,..., q (Graybill, 1983, p. 301). Under case 1 in the example above, where







p = 3, q = 2, o = 1, n = 16, 4K2/n = 64, and m = 12 and if


V=[ .25 -.25
S-.25 .25 '


then

2 64 2 129
varA[f(X)] 2 + -4)r4(.0625 + .0625 + .0625 + .0625) + = 12r + .5.


Recalling that r4 (4 X a)2 = (1 + x1 + X2 + X2)2, we can now determine the average
variance portion of the IMSE as follows:

K f1 11l ft1,j [( + X2 + X2 + x2)2 + .5]dxdx2d3
K 4 varA [ a()]dx = .- l 1 1298"1 1 2 3 + A------
JR f1f1 fl 1, dxedx2dx3

129
= 129 (4.2666667) + .5 = 4.775.
128

The computations above are meant to illustrate that for this 16-run Resolution V
design, the IMSE's can be analytically computed. However, these analytic compu-
tations are not feasible for the other designs being considered so we will numerically
compute the IMSE's using Monte Carlo integration. The results will be displayed
graphically in Section 3.5.

3.4.2 Design 2 and Design 3: 20-run and 24-run Plackett-Burman Designs


Plackett and Burman (1946) demonstrate how to construct Resolution III orthog-
onal arrays involving 4L design runs, L = 2, 3,..., for up to 4L 1 variables. For the
five-variable example, the design matrices for the second and third combined array
experiments were obtained from the first five columns of Plackett-Burman designs
for n=20 and n=24, respectively. These designs are displayed in Figures 3.2 and 3.3
respectively. Note that we did not investigate whether there exists a better choice
of columns from the 20 and 24-run Plackett-Burman designs in their entirety. For




















X1 X2 X3 Z1 z2
1 -1 1 1 -1
1 1 -1 1 1
-1 1 1 -1 1
-1 -1 1 1 -1
1 -1 -1 1 1
1 1 -1 -1 1
1 1 1 -1 -1
1 1 1 1 -1
-1 1 1 1 1
1 -1 1 1 1
-1 1 -1 1 1
1 -1 1 -1 1
-1 1 -1 1 -1
-1 -1 1 -1 1
-1 -1 -1 1 -1
-1 -1 -1 -1 1
1 -1 -1 -1 -1
1 1 -1 -1 -1
-1 1 1 -1 -1
-1 -1 -1 -1 -1


Figure 3.2. 20-run Plackett-Burman Design


















Xl X2 X3 Z3 Z2
1 -1 -1 -1 -1
1 1 -1 -1 -1
1 1 1 -1 -1
1 1 1 1 -1
1 1 1 1 1
-1 1 1 1 1
1 -1 1 1 1
-1 1 -1 1 1
1 -1 1 -1 1
1 1 -1 1 -1
-1 1 1 -1 1
-1 -1 1 1 -1
1 -1 -1 1 1
1 1 -1 -1 1
-1 1 1 -1 -1
-1 -1 1 1 -1
1 -1 -1 1 1
-1 1 -1 -1 1
1 -1 1 -1 -1
-1 1 -1 1 -1
-1 -1 1 -1 1
-1 -1 -1 1 -1
-1 -1 -1 -1 1
-1 -1 -1 -1 -1


Figure 3.3. 24-run Plackett Burman Design







the purposes of this dissertation, a design was chosen only for its ability to estimate

the model. The question of whether a better set of columns exists for estimating the

models of interest (Draper and Lin, 1990) will be left for future research. For these

designs, the P'P matrices are invertible but not diagonal and thus the Cij matri-

ces are not 0, making it impractical to reduce equations (3.7) and (3.9) into simpler

forms.

3.4.3 Design 4: 24-run Orthogonal Main-Effect Plan


Orthogonal main-effect plans represent a class of designs for planning industrial

experiments with mixed levels. Methods for construction of these designs can be

found in many articles including Wang and Wu (1991), Cheng (1989), and Nigam

and Gupta (1985). The 24-run orthogonal array design that we used was obtained

from columns 2, 4, 6, 18, and 20 of the matrix


[B B2 B2]
3B1 B2 -B2


where B = [B1 : B2] is a Hadamard matrix of order twelve (H12) with the first column

of ones removed, B1 is the second column of H12, and B2 represents the remaining

10 columns of H12. The design matrix is illustrated in Figure 3.4.

3.5 Results and Discussion


For each of the four designs and six states of nature corresponding to the example

discussed above, the IMSEs will be calculated using Monte Carlo integration with

points generated over the cuboidal region R, defined in coded variables by R = {x =

(x, x2, x3)': -1 < xi < 1, i = 1,2,3}. Since the designs employ different numbers
of design runs, the IMSEs shall be displayed first unweighted and then weighted


















X1 X2 X3 Z1 Z2
-1 -1 -1 -1 -1
1 1 -1 1 -1
-1 1 1 1 -1
1 -1 -1 -1 -1
-1 1 -1 1 1
1 -1 1 1 -1
-1 1 1 -1 -1
1 -1 1 1 1
-1 -1 1 -1 1
1 1 1 -1 1
-1 -1 -1 1 1
1 1 -1 -1 1
-1 -1 -1 1 1
1 1 -1 -1 1
-1 1 1 -1 1
1 -1 -1 1 1
-1 1 -1 -1 -1
1 -1 1 -1 1
-1 1 1 1 1
1 -1 1 -1 -1
-1 -1 1 1 -1
1 1 1 1 -1
-1 -1 -1 -1 -1
1 1 -1 1 -1


Figure 3.4. 24-Run Orthogonal Array







(multiplied) by the number of runs in order to provide fair comparisons across the

designs.

Figure 3.5 considers the bias contribution to the mean squared error for cases

1-3 discussed above. (Cases 4-6 need not be considered since the bias contribution

is independent of the value of the constant K2 in equation (3.10).) The left hand

column displays the unweighted average bias while the right hand column displays

the weighted average bias. The graphs in Figure 3.5 illustrate the effect of correlation

on the average bias for the example and designs discussed above. By comparing the

vertical axes of the graphs, it can be noted that the average bias increases as the

value of oar = var(zl) increases, although the relative ranking of design performance

is not affected by the value of (11. Note that for the Resolution V design and the

24-run Plackett-Burman design, the average bias is constant over the values of p12

between -1 and 1. Also, when multiplied by the number of runs, the average bias for

the Resolution V design is smallest. With regard to the Plackett-Burman designs,

the 24-run has the smallest average bias when unweighted by the number of runs, but

after weighting, there is little difference between the performance of the 20 and 24-run

designs across the three cases. Interestingly, when weighted, the 20-run is worse than

the 24-run design only for large positive correlation. The orthogonal main-effect plan

had the largest average bias for each case, regardless of weighting. One interesting

feature of the orthogonal main-effect plan however is that the average bias decreases

as the P12 goes from -1 to 1. This suggests the possibility of the existence of a

nonorthogonal design with smaller average bias than an orthogonal design.

Figures 3.6 and 3.7 illustrate the effect of correlation on the average MSE. In

3.6 where the average MSE is unweighted by the respective numbers of design runs,

observe that the 24-run Plackett-Burman design performs best (i.e., has the smallest

IMSEs). Interestingly, in the first and second plots in the left hand column, observe

that the value of aon has an impact on the performance of the Resolution V design










0 = a2 = 0.25


-1.0 -0.5 0.0 0.5 1.0
P 12


a11 =1.0; 022=0.25


1.0 .5 0.0 05 1.0
Pl2


S,,= 2.25; oa =0.25


-1.0 -0.5 0.0 0,5 1.0
P12


,11 = 1.0; 22= 0.25


-1,0 -05 0.0 05 10
P 12


G11=2.25; a22 =0.25
1


-1.0 5 0.0 0,5 LO


-1.0 .5 0.0 05 1O


Legend: RESV .......... PB-20 PB-24 __ OA-24


Figure 3.5. The Effect of Correlation on the Average Bias (Unweighted Weighted)


o = = 0.25


<0










S,,=O 22= 0.25; K = 1152


-10 -05 00 05 1.0



o, =1.0;o = 0.25; K= 128



\ /












10 05 00 0,5 O
P12


T = 2.25; a2= 0.25; K2 = 128




0 05 00 /

N

\ /



., .. -- o' ^ .


\ /












-10 -OS 00 0-5 10






















P 12


a = 2.250 =0.25; c 2= 1152


-10 5 00 0.5 0/
-10 -05 0,0 05 10


Legend: RESV ............ PB-20 ------ PB- 24 OA-24


Figure 3.6. The Effect of Correlation on the Average MSE


o, = Y = 0.25; Kx = 128










an= a22= 0.25; K= 128

-\

\ /
S\ /



\ /
", N / ..,

"',. ___ ,7





-10 0.5 00 0.5 1.0
P 12


ao = 1.0; 22 =0.25; 2= 128

\
-\ /

\/
\ /









.1.0 0.5 00 0.5 1.0
P12




N\

\ /
\ /
/



"'-. -, .. ..- "


-10 -0.5 00 05


1.O0


ao = 22= 0.25; 2= 1152



.\ /




S12






-1.0 -0.5 0,O 0.5 10
P 12


ao =1.0;o5 = 0.25; 12= 1152


0\ /

\ \ /



7





.1,0 -0.5 0.0 0.5 10
P12


oa|=2.25;o22 =0.25; K2 1152


\ /













-1.0 -0,5 0.0 0.5 1.0


Legend: RESV ............ PB-20 PB-24 _- OA-24


Figure 3.7. The Effect of Correlation on the Weighted Average MSE







relative to the 24-run orthogonal array. Also, a comparison of the first plot on each of

the left and right side reveals that the value of K2 has an impact on the performance

of the Resolution V design relative to the 24-run orthogonal array and the 20-run

Plackett-Burman design. Otherwise, the remaining plots show that the values of an

and K2 have little impact on relative design performance. Lastly, note that the average

MSE increases as the value of p12 deviates from zero in either direction. Overall, the

correlation structure of the noise variables has an impact in an absolute sense but

only a small impact in a relative sense. In Figure 3.7, where the average MSEs are

multiplied by the respective numbers of design runs, observe that the Resolution

V design now performs uniformly best while the 24-run orthogonal array performs

uniformly worst. Also, observe that the values of an and K2 no longer have an impact

on relative design performance. The reason for this is that regardless of the value

assumed by K2 = |A12, the variance contribution to the MSE dominates the bias

contribution. Thus, the actual value of JAI seems to have little practical influence on

our recommendations.

In some sense, these results are not very surprising. The superiority of the Res-

olution V design in terms of the mean squared error criterion is enhanced by the

fact that it supports the estimation of the first order C x N interactions in the model.

The latter three designs, however, are strictly main effect plans whose columns were

chosen to yield a nonsingular P'P matrix. The choice of columns allowed the model

to be estimated by the designs, but not as well as the Resolution V design.

3.6 Other Measures of Design Performance


The integrated or average mean squared error criterion combines two aspects of

design performance: average squared bias the average squared distance of f(x) from

the true value of r(x), and average variance the average precision with which -(x)

estimates Tr(a), where average means averaged over the region of interest R. Rather







than considering average design performance, the practitioner may be interested in

the maximum mean squared error of (ax) within the region of interest for a given

design. In this case we seek a design that minimizes the maximum mean squared

error of T(x).

From equation (3.10) the maximum mean squared error can easily be determined:


maxMSE[f (x)] = max{ [bias((B(x))]2 + varAi['(x)]}
XER 2ER

= max{[tr(VE)]2 + 2tr(VE)2
XER

4K2 2r4
+ --tr[(I 0 x,,)VEV(Iq 0 X')] + -- }.
nA n m

For the simulation study and four designs discussed above, the maximum MSE's

were determined by randomly generating points over the cuboidal region defined by

R = {x = (x1, X2, X3)' : -1 < x, < 1, i = 1,2,3}. The resulting maximums are

displayed in Figures 3.8 and 3.9.

In equation (3.10) recall that the variance portion varA[T(x)] was determined

by averaging the variance of f(x) given by equation (3.8) over all possible As with

length K. Hence the maximum mean squared error given above can be thought of as

a maximum-average, i.e., the maximum MSE over the region R after averaging over

the parameter space of A. Instead of averaging var[af(x)] over all possible A's with

length K, consider maximizing var[f(x)] over all possible A's with length K (Noble

and Daniel, 1977):

varM[T(x)] = max var[f(ax)]

2o4
= max {2tr( VE) + 4A'[Iq x]VEV[Iq0 A + --
IAI=n n m










0,= o22= 0.25; K 2= 128

\
- \ /





0-0 \ / L
\ /



S10 0.25; K / 128
0\ 5 ..



.10 -05 00 05 10
P12


al 2.250; o = 0.25; K= 128


0 -05

/










-10 .05 00 0.5 10
P,2


o, = 2.25; o0= 0.25; K2 128

\







\\- /.





-10 -05 00 0-5 1-0


S=0 22= 0.25; K2= 1152

1\ /
\ /
\ /
\ /
\ /


,.. \ / ..*




.10 -05 00 05 10
P 12


a, = 1.0;o2= 0.25; K2= 1152

\ 7/






N- \ / ."."





*1.O -05 00 0.5 10
P12


0O= 2.25;o02 = 0.25; Ic2= 1152


\ /
^\ /

\ / /
\- /

,, ...*- '7
..\





.10 -0. O0 05 t0


SLegend: RESV ............ PB-20 --. PB-24 -_ OA-24

Figure 3.8. The Effect of Correlation on the Maximum-Average MSE









, ,=a22= 0.25; K 2= 1152


\/
\ /
/






,/ 1



S = 1.0 a =2 025 .= 128-
-1.0 05 00 0.5 10








P 12

=2.250;o = 0.25; K= 128
\ /
\ /
\
\ //






-1.0 -0.5 00 0.5 10
P,2

0,,= 2.25;o22=0.25; 12= 128
\ -
S/
\ /

\ /
\ /
\ .




-1.0 -0.5 0.0 0.5 1.0
P12


/
\ /0




P /1






ai 1.0; = 0.25; 2= 1152
-1.0 0.5 0.0 0,5 10








P 12
o,, = 1.0;o=220.25; K2=1152
\ F-
\ /
\ /
N /



J -J



-1O -0.5 0.00 10
P 12

o,,=2.25;022 = 0.25; K 2= 1152



-/

-- I -






P,2


Legend: RESV ............ PB-20 PB-24 OA-24

Figure 3.9. The Effect of Correlation on the Weighted Maximum-Average MSE


o = 22= 0.25; K = 128


I







where em denotes the largest eigenvalue of (Iq 0 )a)VEV(Iq x'). This yields a
maximum-maximum mean squared error criterion:


maxMSE[4 (x)] = max{[bias( (x))]2 + varM[(xa)]}
XER 2XR

24
= max{[tr(VE)]2 + 2tr(VE)2 + 42em +-- }.
n-rn


For the simulation study and four designs discussed above, the maximum MSE's

were determined by randomly generating points over the cuboidal region defined by

R = {x = (Xl,X2,X3)' : -1 < xi < 1,i = 1,2,3}. The resulting maximums are

displayed in Figures 3.10 and 3.11.

Figures 3.8 3.11 illustrate that when weighted by the number of design runs,

the Resolution V design provides the lowest values of the maximum-average and

maximum-maximum mean squared errors. Interestingly, when unweighted by the

number of design runs, the Resolution V design still provides the lowest values of

the maximum-average and maximum-maximum mean squared errors except when

P12 = cov(zl, z2) = 0. When the noise variables are independent the 24-run Plackett-
Burman design provides the lowest values. Overall, the Resolution V design per-

forms best under these maximum-average and maximum-maximum criteria whether

weighted or unweighted by the number of design runs. This differs from the conclu-

sions obtained with regard to the integrated mean squared error criterion. In that

case, the Resolution V design was recommended under weighting while the 24-run

Plackett-Burman design was the clear winner when unweighted by the number of de-

sign runs. Lastly, a comparison of the left and right hand columns of plots in Figures

3.8 3.11 indicates that once again the value of r2 (128 or 1152) has little impact on

the relative performance of the four designs.





-1.0 -0.5 0.0 0.5 10
P12


al = 2.25; a22 = 0.25; 2= 128



N /
/

/


S. *\--. ..*"


-10 05 0.0 0.5


10


11 = 22= 0.25; K'= 128
-\
\ /
\ /

\ /
\ / ..;







-10 -0.5 00 0.5 1.0
P 12


a0, = 1.0; a22 = 0.25; K2= 128


/


; //
\ /


"-';.. \ ..' ^
~~\ ':-- / **'/ >


Legend: RESV .......... PB-20 ----- PB-24 _. OA-24

Figure 3.10. The Effect of Correlation on the Maximum-Maximum MSE


0r11= 0 22-= 0.25; K 2= 1152








S = 1.0;22 0.25; ic 2= 1152
\
\ /

\ /
\/\ /



., \ / .. "




-10 -05 00 05 10




S\/ /

\ .' /











-1.0 -0.5 00 0,5 1.0
P.2


Oau = 2.25; o= 0.25; c2 = 1152




\ -//

\- /
\.

,





-1,0 -0,. 0.0 0.5 10











a, = o~= 0.25; = 128


\/











S12





11 1.0; 22 0.25; 2 128
-10 -05 0-0 05 10
P 12


O11 =1.O;o2z = 0.25; Kc 128

\/
\ /
\ /



\ /




i i i

-1 0 -0. 00 05 10
P12


a11=2.25;o=o0.25; K= 128


\ /
N /



0'-05 00 1) 10-



-1 0 .0.5 0.0 0.5 1.0


S11= 22= 0.25; K2= 1152
















-10 _0-5 00 0.5 10
\ /












S225 0.25 1152
00 0
/



S. \ I ,/




-10 -0,5 00 05 10
P,2


ao =21.0;o =0.25; K2 1152

X 7I
\ /
/




\. /-






- 10 .05 00 0.5 0


r 12 r 12

Legend: RESV ............ PB-20 -----. PB-24 -. OA-24


Figure 3.11. The Effect of Correlation on the Weighted Maximum-Maximum MSE







3.7 Practical Application


In this section we will illustrate our methodology for an example discussed in Box and

Jones (1992). They consider a setting in the food manufacturing industry in which

an experimenter is seeking the best recipe for a cake mix to be sold in a box at the

supermarket. For this example, we suppose that there are four control factors, flour

(F), shortening (S), sugar (G), and egg powder (E). Next suppose the experimenter

anticipates that a consumer's oven temperature may be biased up or down and that

a consumer may overcook or undercook the cake. Thus, there are two noise variables,

time (t) and temperature (T), and it is of interest to determine a recipe, i.e., settings

of F, S, G, and E that yield a good-tasting cake even if the consumer strays somewhat

from the baking time and temperature recommendations on the box.

For this combined array experiment having four control and two noise variables,

assume we are interested in fitting the model


y = 3o N+ O13iX + /2X2 + #3X3 + /4X4 + /l1Zl 7222

+ Alxlzli + A21x2z1 + A31X3Z1 + A41x4z1

+ A12XlZ2 + A22X2Z2 + A32X3Z2 + A42X4Z2 + 6.

Let us consider six reasonable experimental designs that can be used to fit this model:

Design 1: 32-run Resolution VI Fractional Factorial Design (26-1),

Design 2: 32-run Taguchi Crossed Array Design (24-1 22),

Design 3: 22-run D-optimal design recommended by Shoemaker et al. (1991)

(this design was intended to be a smaller competitor to Taguchi's 32-run crossed

array and was generated using RS/DISCOVER to fit a model that also included

the interactions X1X2, x1x3, x1x4, and zlz2),







Design 4: 22-run D-optimal Design generated by SAS's OPTEX procedure,

assuming Shoemaker's model,

Design 5: 22-run D-optimal design generated by SAS's OPTEX procedure,

assuming our model, and

Design 6: 16-run D-optimal design generated by SAS's OPTEX procedure,

assuming our model.

In order to compare these designs utilizing the integrated mean square error cri-

terion developed in Section 3.3, the experimenter must specify V and K'2 in equation

(3.10). For this example, let us assume that the recommended baking time and tem-

perature are 40 minutes and 325F, respectively. We anticipate that the consumer

may bake the cake for 35 to 45 minutes at an oven temperature between 300F and

3500F. Hence we may estimate all = var(zi) using (range/4)2 = (2.5 min)2 and

a22 = var(z2) using (range/4)2 = (12.5F)2. Assuming each variable is studied at

two levels, denoted by -1 and +1 in coded design units, then this corresponds to

taking ano = 22 = (.5)2 = .25. Also, it would seem reasonable to assume that in

this example the noise variables are negatively correlated. This is because we be-

lieve that a consumer is likely to know for instance that his or her oven temperature

is biased upward and will try to compensate by cooking the cake for less than the

recommended amount of time. Thus we will consider P12 equal to -0.5, -0.25, and

0 which respectively indicate strong negative correlation, weak negative correlation,

and independence between cooking time and temperature. To specify K2, we recall

from the simulation study that the results were robust to the actual choice of JAI, so

we will assume that the magnitude of each element of A is four times the magnitude

of o,. If we take a, = 1, then 02 Al2 = nA(16ac) = (10)(16) = 160.

The IMSEs for the competing designs were calculated using Monte Carlo integra-

tion with points generated over the cuboidal region R, defined in coded variables by







R = {a = (xi,x2,X3,X4)' : -1 < xi < 1, i = 1,2,3,4}. The results are displayed in

Figure 3.12. Since the designs have different numbers of runs, we display the graphs

first unweighted and then multiplied by the number of runs to allow for fair com-

parisons. The Taguchi design is not present in the graphs since its P'P matrix is

equivalent to the P'P matrix from the fractional factorial design. Thus, the designs

have equivalent C matrices and so the MSEs calculated from equation (3.10) would

be equivalent apart from simulation error.

From Figure 3.12 we note that the 32-run resolution VI fractional factorial design

has the smallest bias and MSE in both the weighted and unweighted comparisons.

Interestingly, Shoemaker's D-optimal design generated using RS/DISCOVER has sig-

nificantly larger bias and mean squared error than the 22-run D-optimal design we

generated using the OPTEX procedure in SAS. Hence, it seems that the algorithm

for generating D-optimal designs can have some impact on our criterion. Another

interesting feature is that the 22-run D-optimal design that was generated to fit our

smaller model performed worse than the two 22-run designs used to fit Shoemaker's

larger model. Lastly, note that the 16-run D-optimal design, the smallest design

that is allowed by our criterion since n m must be greater than zero, performed

significantly worse than the other four designs, regardless of weighting. After remov-

ing the 16-run D-optimal design from consideration, Figure 3.13 provides a better

comparative display of the remaining four designs.

3.8 Conclusions


The combined array is an experimental design used for finding conditions that

achieve a target condition for a response of interest while also minimizing squared

error loss. Under the assumption of an additive error structure, the appropriate

squared error loss involves the estimated process variance.



















O = 22= 0.25


-. ------ -





















-0, -04 -03 402 -01 00




5 -04 -,3 -02 01 010
P12


In = o22 = 0.25
















S -04 3 -02 -01 00
P 12


a I= 22= 0.25
















-05 -04 0 02 4O1 00
P 12


Legend: DOPT16 ............ SHOE22 --- DOPT(SHOE22)

DOPT22 RES-VI


Figure 3.12. The Effect of Correlation on the Average Bias and Average MSE (Un-
weighted Weighted)
























0 =11 C 22= 0.25


-- ------------- ....

- -- -


-0. -0 -03 0.2 -01 0,0

P12



011= 022= 0.25

--. --.-


-0.5 -0.4 03 -02 01 0.0
P12


0.=0 22= 0.25

-- .. ... .. . .


-05 04 03 -02 0.1 00

P 12


n= g 22= = 0.25


N. -







N "-. "- 2... ..
^ '''--------'-- ^


-05 -04 03 02 -01


00


Legend: ........... SHOE22 -----. DOPT(SHOE22) DOPT22 RES-VI


Figure 3.13. The Effect of Correlation on the Average Bias and Average MSE (Un-
weighted Weighted)


5,

~







Our research establishes that the combined array yields a biased estimate of the

loss. This bias is maximized at the perimeter of the region of interest and can be

relatively large if the design is highly fractionated, the noise factors are highly variable

in the process, or the inherent variability is large.

In terms of the mean squared error of the estimated loss, our research indicates

that the variance of the estimated loss tends to dominate the bias. We have shown

that the nature of the correlation structure can have a significant impact on the mean

squared error. However, the relative comparison of competitor designs through our

mean squared error criterion was impacted very little by the choice of V, the variance-

covariance matrix of the noise factors, or n, the length of the vector of coefficients.

On the whole, when the number of design runs is taken into consideration, Resolution

V and higher fractions of 2k factorial designs appear to produce estimated losses with

the smallest mean squared error.












CHAPTER 4
ANALYSIS OF THE COMBINED ARRAY
WHEN THE ASSUMED MODEL IS UNDERSPECIFIED


4.1 Introduction


In Chapter 3, the properties of the estimated process variance and the resulting

conclusions regarding various combined arrays were derived under the assumption

that the presumed model for the response was correct. These derivations represented

the first formal attempt at considering the properties of the combined array. Com-

bined arrays are currently being used in practice but little attention has been devoted

to the properties of the estimator of the process variance. The importance of consid-

ering the properties of the process variance estimator lies in the fact that the analyst

considers predicted values of the process variance to determine the robust control

factor settings. We showed in Chapter 3 that the estimator of the process variance

is biased even though the presumed model was correct. This led to the derivation of

a mean squared error criterion which was used to determine which combined arrays

among a given set produced process variance estimators having the smallest mean

squared error for a particular experimental situation.

It is widely recognized that the fitted model serves only as an approximation

of the true model and the practitioner prefers the fitted model to be as simple as

possible. Hence the fitted model is often times an underrepresentation of the true

response surface. Thus, it is of interest and importance to determine the impact of

underspecifying the model on the properties of the estimated process variance. Using

the foundational notation and methodology that have been set up for the correct







model case, in this chapter we will extend our results to the case where the model is

presumed to be underspecified. Our goal once again is to derive the mean squared

error of the estimated process variance in order that an integrated mean squared error

criterion can be used to evaluate the performance of combined arrays. Ultimately, we

want to specify which combined arrays yield estimates of the process variance that

are least affected when the underlying model may have been underspecified.

We will begin by formulating the problem of model misspecification in its most

general setting, and in the next chapter discuss special cases of interest within the

details of various examples.

4.2 Properties of the Estimated Process Variance


Suppose an experimenter fits the following model to a quality characteristic of

interest:

y = 00 + x'i31 + z'y71 + x 'lAlli + e (4.1)

where

/#o is the intercept,

x, is a pi x 1 vector representing the appropriate polynomial expansion of the

control variables of order dc, for the assumed model,

$1 is the pi x 1 vector of unknown parameters associated with the control

factors,

zl is a ql x 1 vector representing the appropriate polynomial expansion of the

noise variables of order ds, for the assumed model,

-y is the q, x 1 vector of unknown parameters associated with the noise factors,

Aln is the pi x ql matrix of unknown parameters associated with the control-

by-noise interactions (some elements of All may be zero), and







e is the random error.

However, let 7r be the true response represented by

77= 03o+x~ +x2I32+z~7l1+z2'yz2+'Anl+xlA12z2+a2A21z1+zx2A22z2+E (4.2)

where /o, a1, i1, z3 ~y1, A11, and e are as defined above and

X2 is a p2 x 1 vector representing the appropriate polynomial expansion of the
control variables of order dC2,

f32 is the associated P2 x 1 vector of unknown parameters,

z2 is a q2 x 1 vector representing the appropriate polynomial expansion of the

noise variables of order dn2,

Y2 is the associated q2 x 1 vector of unknown parameters,

A12 is the matrix of unknown coefficients of the control-by-noise interactions

between the control variables in xa and noise variables in z2,

A21 is the matrix of unknown coefficients of the control-by-noise interactions

between the control variables in x2 and noise variables in z1, and

A22 is the matrix of unknown coefficients of the control-by-noise interactions

between the control variables in a2 and noise variables in z2.

Note that the results that follow in this chapter do not depend on the orders of the

polynomial expansions of the control and noise variables in the fitted and true models.

However, the smallest model underlying a combined array experiment must contain

first order terms in the control (C) and noise (N) variables as well as first order

control-by-noise (CxN) interactions. Only by exploiting the CxN interactions can

the experimenter determine the effect that changing the levels of the noise variables

has on the control variables. Hence, some but not all of the elements of Aln may







be assumed a priori to be zero in the fitted model, while some or all of the elements

of A12, A21 and A22 may be taken to be zero in the true model. In Chapter 5 we

will take d,1 = dc, = 1 and d, = d,, = 2 in order to investigate the impact of

underspecification in the control variables alone, noise variables alone, control-by-

noise interactions alone, and some combinations of the three.

For models (4.1) and (4.2) we shall assume that the es are independent and iden-

tically distributed with mean zero and variance or. For the purposes of a particular

planned experiment, the noise factors shall be treated as having fixed effects. How-

ever, in the process we shall assume zl follows a multivariate normal distribution with

E(zi) = p, and var(zi) = V1. Likewise we assume z2 follows a multivariate normal

distribution with E(z2) = 12 and var(z2) = V2. We shall assume aui = var(zi)

and oij = cov(zi, Zj) are known and depending on the order of the vectors zi and

z2 additional elements of Vi and V2 may be derived through differentiation of the

moment generating function of the multivariate normal distribution. These compu-

tations are shown for a special case in Appendix D. Furthermore, the error terms and

noise variables shall be assumed to be independent.

To develop the prediction properties of the estimated process variance, it will be

helpful to establish a set of notation which follows:

A'/Ia =[71,,A'1]

S= (1,'1)

A11 = vecAlna



1=1 (00,"1,Al)

A12a = A12, A2]

A12 = vecAi2a

A21 = vecA21







A22 = vecA22

'2 = ( 2 12 ) 211 \22

02(Xi) = 2 + A'121 = A'2aXia = (Iq2 0 Xja)A12

A212 = ( 0 2)21

A2 = (q2 0 X)A22

n = total number of design runs

mi = total number of unknown coeficients in model (4.1) = (1 + pi)(l + qi)

M2 = total number of unknown coeficients in model (4.2) = (1 + p1 + P2)(1 + q, + q2)-

Using the above notation, models (4.1) and (4.2) can be expressed in the more

general forms given by:
y = P101 + e (4.3)

and
ri = P PI 2 + -2 2 + e (4.4)

where

P1 is the n x mi (mi = 1 + pi + ql(pi + 1)) model matrix whose elements are

known functions of the settings of k input variables determined by an n x k

design matrix, and

P2 is the n x (m2 m1) (m2 ml = p2 + q2(P2 + 1)) matrix associated with

terms not present in the fitted model but present in the true model.

Assume the appropriate combined array experiment has been run to estimate

the linear model given by (4.3). Using ordinary least squares to estimate the model

coefficients, the resulting estimate of the coefficient vector is i1 = (PP1i)-P1y

Now if the true model is given by (4.4) and we assume the noise variables to be fixed

in the experiment, then E[y] = P101 + P202. Thus the prediction properties of 4i






are as follows:

E[] = E[(P'P)-1 P'y]

= (P'Pi)-'PIE[y]

(PIP)-lP(PI'1 + P212)

= (PlPI)-1PlPI-,1 + (PIP1) PIP2b2

= 01 + A 2 (4.5)

where A = (P'P1) -P'P2 is called the alias matrix and

var[ ] = P)-1

Closer consideration of equation (4.5) will enable us to derive the prediction properties
of A11:



E[I]= E(qI) = (3 A A12
E(AI) A11 A 2 1
\ 22 /

Thus, it follows that

E[An1] = All + A*02 (4.6)

where A* is formed by deleting the first p, + 1 rows of the alias matrix A. Now, since

k1(Xi) = (Iqi 0 x31)A1, it follows from equation (4.6) that

E[ 1(xi)] = (Iq 'ia)(A11 + A*02)

= 41(z1) + (Iq 0 .ia)A*I2. (4.7)

We are now in a position to derive the prediction properties of the estimated
process variance which we denote by ,(xl). If the fitted model is given by (4.1) and
the noise variables are truly random in the process, then analogous to equation (3.2)







the estimated process variance may be expressed as

(Mxi) = var[y] = 4~(xl)Vpk(xi) + &o.

If the true model is given by (4.2), the derivation of the expected value of Ty(xa)
follows from (4.7)

E[,ry(ax)] = tr(ViSi) + E[ j(xl)]'ViE[1l(xl)] + o'

= tr(vi i) + [i('i) + (Iq 0 ,o)A* '']V1[1(,1i) + (Ij xi)A*l2] + 02

= tr(Vix) + 0Y1(Xi)V1i1(ai) + 204((1)V1(Iq, 0 a'j)A*'

+ 2A*(IJ, 0 Xla)Vl(Iql, la)A*02 + 0 (4.8)

where E, = var[ 1(xi)]. Also, from (4.2) we can derive the true process variance,
which we denote by r,((x), x = (x'1,x')', as follows:

-,(x) = var[y] = var[(i' + -xA11 + 'A21)z1 + (72 + ,'A12 + x'A22)z2 + E]

= (7' + xAnll + x'2A21)Vi(7- + BAllA + '2A21)'

+ (7- + 'A12 + x'A22)V2(Y72 + A12 + x'A22)'

+ 2(7'i + An11 + x2A21)V12(72 + 1xA12 + x'A22)' + o

= 4'(xi)V11(xi) + 20'(x1)VIAa2 + X1A21VA'21 2

+ '2(al)V202(al) + 24'(aX)V2A2x2 + 2xA22V2A'222

+ 2(0((xi) + x' A21)V12( 2(X) + A'222) + U2 (4.9)

where V12 = cov(z1,z2) and the terms 0'(xa) and '2(x1) are explained in the list
of notation given on pages 62 and 63.







4.3 Expressions for the Bias and Variance of the Estimated Process Variance


From (4.8) and (4.9) note that the estimated process variance is biased. In Ap-

pendix B it is shown that the bias can be conveniently expressed as

bias[fy(axi)] = E[f,(xi)] r,(x)

= tr(V1 E) + a'Ba (4.10)

where a' = (A,, ip) and the matrix B is is given by equation (B.3).

In order to develop a mean squared error criterion for estimating the loss incurred

by using ,(axi) to estimate r,(ax), we must also derive the variance of fy(xi). If in

model (4.4) we assume that e is distributed as multivariate normal with mean 0 and

variance oIr, then in Appendix C it is shown that the variance of the estimated loss

can be conveniently expressed as


var[y,(xi)] = 2tr(ViEl)2 + 4a'Wa + 2a (4.11)
n m1

where a'= (A', i') and the matrix W given by equation (C.2).

Both the bias and the variance of Ty(x1) given by equations (4.10) and (4.11)

respectively depend on the unknown parameter vector a. We proceed here as in

Chapter 3. Let us assume that the length of the vector a is known where K = lal.

Define U to be the surface of the unit sphere and ds, to be a differential on this

surface. Then by averaging the bias and variance over all possible as with length K

it is possible to eliminate the dependence on a. The averaged bias and variance are

given by

2 KK2
biasA[y(xl)] = K2K bias[f,(xai)]ds, = tr(ViSi) + -trB (4.12)
JU n,

and

varA[y(Xi)] = 2K J var[fy(xi)]dsa
JU







4K2 2U4
= 2tr(ViEi)2 + -trW + (4.13)
no n mi

where K = [fu dsa]-', and n, is the number of elements in a not assumed a priori to

be zero.

Lastly, from equations (4.12) and (4.13) we obtain the following equation for the

average mean squared error of Ty(xi) which can be readily calculated at any point

X1:

K2
MSEA['y(xl)] = [tr(ViE) + -trB]2

4K2 294
+ 2tr(ViEi)2 + -tW + -- (4.14)
n, n mi

4.4 The IMSE Criterion and Other Measures of Design Performance


As discussed in Chapter 3, a reasonable method for comparing designs is to eval-

uate the integrated mean squared error (IMSE) of Ty(xl) for each design. The IMSE

is given by

IMSE[fy(a)] = K MSEA[?y(x1)]dx


where R is the appropriate region of interest, K = [fRdx~]-', and MSEA[y1(xa)] is

given by equation (4.14).

It may also be of interest to compare designs by evaluating the maximum mean

squared error of N(xE) which is given by


max MSEA[f'y(Xl)] = max{[biasA(,y(Xl))]2 + varA[fy(X1)]}.
XiER X1ER

Since MSEA[y1(xl)] was determined only after averaging bias[f',(xi)] and var[fy(xi)]

over all possible as with length K, this maximum mean squared error can be thought

of as a maximum-average. An alternative maximum mean squared error criterion can







be determined by first maximizing the bias and variance of fy(ax) over all possible

as with length K as follows:


biasM[y(xl)] = max[tr(ViEi) + a'Ba]


= tr(ViEl) + K2eM(B)


where eM(B) denotes the largest eigenvalue of B, and

24
varM[Ty(i)] = max[2tr(ViEi)2 + 4a'Wa + --
laI=K n mi


= 2tr(Vl~l)2 + 4K2eM(W) + -2
n m1

where eM(W) denotes the largest eigenvalue of W. This yields a maximum-maximum

mean squared error criterion given by


max MSEM['y(xl)] = max{[biasM(jy(xl)]2 + varM[y (a )]}.
X1ER XiER

We are now ready to extend the case study and practical example that were discussed

in Chapter 3 to the situation where the presumed model is underspecified.












CHAPTER 5
SPECIAL CASES OF MODEL UNDERSPECIFICATION


5.1 Introduction


In this chapter we consider special cases of model misspecification. For all cases

we will begin by assuming that the fitted model is first order in the control (C) and

noise (N) variables and contains first order control-by-noise (C x N) interactions. In

other words, in model (4.1) we take x1 = (xi,..., ,Xp)' to be a vector of order d,, = 1

and z1 = (1, ,. .. zq)' to be a vector of order d, = 1. We will consider eight different

"true" models containing the same terms in the fitted model plus additional higher

order terms not present in the fitted model. Essentially, we assume the fitted model

is underspecified. The additional higher order terms that will be incorporated into

the true model are as follows:


Case 1 (C2): The true model contains pure quadratic terms in the control

factors and first order control-by-control interactions.

Case 2 (N2): The true model contains pure quadratic terms in the noise factors

and first order noise-by-noise interactions.

Case 3 (C2 x N): The true model contains interactions between noise factors

and second order control effects.

Case 4 (C x N2): The true model contains interactions between control factors

and second order noise effects.







Case 5 (C2, N2): The true model contains pure quadratic terms in the control

factors, first order control-by-control interactions, pure quadratic terms in the

noise factors, and first order noise-by-noise interactions.

Case 6 (C2, C2 x N): The true model contains pure quadratic terms in the con-

trol factors, first order control-by-control interactions, and interactions between

noise factors and second order control effects.

Case 7 (N2, C x N2): The true model contains pure quadratic terms in the noise

factors, first order noise-by-noise interactions, and interactions between control

factors and second order noise effects.

Case 8 (C2, N2, C2 x N, C x N2): The true model contains pure quadratic terms

in the control factors, first order control-by-control interactions, pure quadratic

terms in the noise factors, first order noise-by-noise interactions, interactions

between noise factors and second order control effects, and interactions between

control factors and second order noise effects.


These eight cases cover a variety of experimental situations that the practitioner may

be concerned about. Under differing circumstances the experimenter may fear the

presence of higher order terms in the control variables alone (Case 1), noise variables

alone (Case 2), control-by-noise interactions alone (Case 3 and Case 4), or some

combination of higher order terms in both the control and noise variables (Cases

5-8).

5.2 Reduced Expressions for the Bias of the Estimated Loss


For each special case, the general expression for the bias of the estimated loss that

was derived in Appendix B can be reduced into a simpler form. The simplification is

a direct consequence of setting appropriate terms equal to zero in the bias expression.







The elements of the B matrix in equation (B.3) depend on Q1, Q2,-..., Q1 which
are matrices of the quadratic and bilinear forms described in Appendix B. For the
fitted and true models that we are considering, note that Q8, Q9, Q10, and Q11 are
null matrices. This first simplification is a consequence of the fact that the fitted
model contains only first order noise variable terms and the true model contains at
most second order noise variable terms. In this case V12 = cov(zi, z2) = 0 under the
assumption that z = (z1, z2)' follows a multivariate normal distribution with mean
JA and variance-covariance matrix

v=(VIl V12)
V'12 V2 "

Since V12 = 0, clearly Q8, Q9, Q10, and Q11 are null matrices.

5.2.1 Case 1 Reductions


In case 1, the true model contains C, N, C x N, and C2 terms. Thus, in equation
(4.2) we set 'Y2, A12, A21, and A22 equal to zero. This model is given by

= /#0+ X'131 + x'232 + z1Y1 + x'1Alz + (5.1)

Under these conditions, any terms involving M1, M2, and M3 would be eliminated,
and so the matrices B1 and B2 in Appendix B reduce to B1 = Q1 and B2 = Q2.

5.2.2 Case 2 Reductions


In case 2, the true model contains C, N, C x N, and N2 terms. Thus, in equation
(4.2) we set 132, A12, A21, and A22 equal to zero. This model is given by

] = 00 + X111 + z'1iY + z1y2 + 'AnAzi + e. (5.2)

Under these conditions, any terms involving M2 and M3 would be eliminated, and
so the matrices B1 and B2 in Appendix B reduce to B1 = Q1 and B2 = Q2 Q5.







5.2.3 Case 3 Reductions


In case 3, the true model contains C, N, C x N, and C2 x N terms. Thus, in
equation (4.2) we set /2, Y2, A12, and A22 equal to zero. This model is given by

S= 0o+ x31 + z'iy1 + xAlzzi + x'2A21z + e. (5.3)

Under these conditions, any terms involving M1 and M3 would be eliminated, and so
the matrices B1 and B2 in Appendix B reduce to B1 = Q1 Q3 and B2 = Q2 Q4.

5.2.4 Case 4 Reductions


In case 4, the true model contains C, N, C x N, and C x N2 terms. Thus, in
equation (4.2) we set /2, 'Y2, A21, and A22 equal to zero. This model is given by

-7 = o + a'/31 + z'-y + a'A-nzi + x'A12z2 + e. (5.4)

Under these conditions, any terms involving M2 and M3 would be eliminated, and
so the matrices B1 and B2 in Appendix B reduce to B1 = Q1 and B2 = Q2 Q5.

5.2.5 Case 5 Reductions


In case 5, the true model contains C, N, C x N, C2, and N2 terms. Thus, in
equation (4.2) we set A12, A21, and A22 equal to zero. This model is given by


q1 = /o + x/31 + X'202 + z'721 + z2'72 + x'lAliz + e. (5.5)

Under these conditions, any terms involving M2 and M3 would be eliminated, and
so the matrices B1 and B2 in Appendix B reduce to B1 = Q1 and B2 = Q2 Q5.







5.2.6 Case 6 Reductions


In case 6, the true model contains C, N, C x N, C2, and C2 x N terms. Thus, in
equation (4.2) we set 72, A12, and A22 equal to zero. This model is given by

=0 + X'131 + x2'32 + z[1 + xaAniz + x2A21 z + (5.6)

Under these conditions, any terms involving M1 and M3 would be eliminated, and so
the matrices BR and B2 in Appendix B reduce to B1 = Q1 Q3 and B2 = Q2 Q4.

5.2.7 Case 7 Reductions


In case 7, the true model contains C, N, C x N, N2, and C x N2 terms. Thus,
in equation (4.2) we set /32, A21, and A22 equal to zero. This model is given by

0= /0 + X'11 + z' Y1 + z2'2 + zAnzllz + zxA12 z+ (5.7)

Under these conditions, any terms involving M2 and M3 would be eliminated, and
so the matrices B1 and B2 in Appendix B reduce to B1 = Q1 and B2 = Q2 Q5.

5.2.8 Case 8 Reductions


In case 8, the true model contains C, N, C x N, C2,N2, C2 x N, and C x N2
terms. Thus, in equation (4.2) we set only A22 equal to zero. This model is given by

i = o/ + 1301 + X'202 + Z7 11 + z"'y2xAzi + x'2A21Z1 + a 'A12z2 + e. (5.8)

Under these conditions, any terms involving M3 would be eliminated, and so the
matrices B1 and B2 in Appendix B reduce to BI = Q1 Q3 and B2 = Q2 Q4- Q5.

5.3 Example


In Chapter 3, we noted that even though the fitted model was presumed to be
adequate, the process variance estimator was biased. The development of a mean







squared error criterion for the estimated process variance and subsequent simulation
study enabled us to compare a set of designs. In this section, where the fitted model

is presumed to be underspecified, we will conduct another series of simulations to

investigate the additional impact of model misspecification on the bias and variance

of (xzi). For each of the eight special cases we will again consider a robust parameter

design experiment consisting of three control variables and two noise variables, each

being studied at two levels. For this example, the fitted model given by (4.1), in

expanded form is:

y = 3o + /1x1 + #2x2 + 03X3 + 71zi + 7y22 + A11xzi + A21x2zI + A31 x3z1

+ A12X1Z2 + A22X2Z2 + A32X3Z2 + E.

We will compare the same four designs utilized in Chapter 3 to fit this model: a 16-run

Resolution V fractional factorial design, 20 and 24-run Plackett-Burman designs, and
a 24-run orthogonal main effect plan. We will compare these designs by evaluating the

integrated mean squared error of B(xi), and the maximum-average, and maximum-

maximum criteria that were derived in Section 4.4. These criteria will be calculated

using Monte Carlo techniques with points generated over the cuboidal region defined

by R = {x= (x,x2,x3)' :-1 < xi < 1,i= 1,2,3}.

Recall from equation (4.14) that computation of the mean squared error of 'y(x()

at any point xla depends on V1, the assumed known variance-covariance matrix for

the noise variables in the fitted model, and the value of K2/n,. For this example with

two noise variables in the fitted model, let j1i, a12, and a22 denote the elements of

the V1. Using the same argument from Chapter 3 based on specifying the range of

the noise variables, we will consider three cases:

1. anl = a22 = .25 (R1 = Ri, R2 = );


2. all = 1; a22 = .25 A (t1 = Ri/2, R2 = R2);






3. all = 2.25; 022 = .25 # (iR = Ri/3, R2 = R2);

for which the corresponding VI matrices are given by

1P 1 1 P12 I 1 2 P12 and V, = 3 [ 3 P12
V=4 P12 1 j' 12 P12 1/2 and 4 P12 1/3


To investigate the effect of the correlation structure of the noise factors on the bias
and mean squared error, we will consider a sequence of values of P12 between -1 and
1. Lastly, note that in special cases 2, 4, 5, 7, and 8, where the vector z2 of second
order noise effects is present in the true model, we must determine V2, the variance-
covariance matrix of these second order effects. The general form of V2 is derived in
Appendix D.

5.4 Results and Discussion


For each special case of the underspecified model we made use of the simplifications
discussed above and used SAS/IML to compute the average bias and average mean
squared error of T,(xai) for the four designs under consideration. Eight special cases,
two values of K2, and three Vi matrices produced a total of forty eight sets of values.
As in Chapter 3, we will display the biases and mean squared errors graphically, both
unweighted and weighted (multiplied) by the respective numbers of design runs, in
order to make conclusions and recommendations regarding the four designs. The
graphical displays will illustrate the impact of: (1) the variances, all and U22, of the
noise factors; (2) the correlation, P12, between the noise factors; and (3) the length of
the vector a (2 = a 2). These plots are given by Figures 5.1 through 5.32. When
observing the graphical displays, bear in mind the three-fold intent of this simulation
study:

1. To investigate whether the average bias and average mean squared error are
"robust" to the values of '11, 0"22, P12, and K2.







2. To investigate the impact of model underspecification on the bias and mean

squared error of (axi).

3. To determine which design performs best within each case and across cases of

the underspecified model. By "performs best", we mean that the design has

the smallest bias and/or mean squared error of fy(xi).

Figures 5.1 and 5.2 display the average bias for case 1 in which the true model

contains additional pure quadratic and interaction terms in the control variables only.

In Figure 5.1, where the average bias is unweighted by the number of design runs,

observe that the 16-run Resolution V fractional factorial design (RESV16) performs

best by a small margin over the 24-run orthogonal array (OA24). It is important to

note that the difference in the average bias associated with the RESV16 and OA24

designs is significant based on a 95% confidence interval computed using simulated

standard errors. The scale on the vertical axis causes the difference to appear non-

significant. The 24-run Plackett-Burman design (PB24) is third best followed by

the 20-run Plackett-Burman design (PB20) in last place. One interesting feature to

notice is that the average bias for the PB20 decreases as p12 ranges from -1 to 1,

while the average bias is nearly constant for the other three designs over the values of

P12. Comparing the left and right columns of plots we note that apart from a change

in scale the relative performance of the four designs does not depend on the value

of K2. In Figure 5.2 where the average bias is multiplied by the number of design

runs, we see that necessarily the RESV16 still performs best with the OA24 coming

in second. Interestingly, the PB20 design performs worst when the noise variables

are negatively correlated while the PB24 performs worst when the noise variables are

positively correlated.

Figures 5.3 and 5.4 display the average mean squared error for case 1. Interestingly,

in Figure 5.3, where the average MSE is unweighted by the number of design runs,

we see that the RESV16 performs best with one exception. When K2 = 128, acr =









o 1 = a = 0.25; K = 1152


10 -0.5 0.0 0.5 10
S12


a = 1.0; a22 = 0.25; K1= 128




- - -


1.0 05 00 05 10
P12


a0I = 2.25; a22= 0.25; i =_ 128


............--*....


1.0 -0.5 0.0 0.5 10
P 12


0 = 1.0; 22= 0.25; Ki 1152


............------- .......


Ii
-1o0 405 0.0 0-5 10
P 12


o ,= 2.25;Y22 =0.25; K = 1152

... . -


-10 0.5 00
P 12


05 1.0


.10 0.5 00 05 10


Legend: RESV ............ PB-20 ------ PB-24 OA-24


Figure 5.1. Case 1 (C2) The Effect of Correlation on the Average Bias


o, = y,2 = 0.25; K =_ 128









0 = 0=2 0.25; = 128

. .- :-- _- -


10 .0- 00 0.5 10
S12


0,, = 1.0; o22 = 0.25; )= 128


- -


10 05 00 05 10
P 12


0 I= 2.25;022= 0.25; x = 128


.. . . .


S= ao22= 0.25; ic = 1152














.1.0 -05 00 0.5 10
P 12


o1, = 1.0; 22= 0.25; i= 1152


-- T -.- -


10 -05 00


05 10


P 2


0a1,=2.25;022 = 0.25; ie= 1152


10 0.5 0,0
P12


0,5 10


1,0 -05 00 05 10


OA-24


Figure 5.2. Case 1 (C2) The Effect of Correlation on the Weighted Average Bias


Legend: RES ............ PB-20 --. --. PB-24











o = 22= 0.25; = 128


" .......... .








-LO -0,5 0.0 0.5 10
P 12


a = 1.0;22 = 0.25; i = 128


.10 .05 00 05 10
P 12



1 = 2.25; a = 0.25; ? = 128





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











I- -. o0. 0.5 1L


0, = 022= 0.25; K 2 1152


-10 -0.5 00 05 10
P 12


ao1 = 1.0; 2= 0.25; i = 1152


O -05 0.0 05 10

P 12


a, = 2.25;oa = 0.25; e = 1152


.i -n5 on 05s O


Legend: RESV ............ PB-20 ----- PB-24 OA-24


Figure 5.3. Case 1 (C2) The Effect of Correlation on the Average MSE


i0


- -- -- --- --











0 = o,, = 0.25; c 1152


10 -0.5 00 05 10
P 2


0,, =1.0; 0 =0.25; = 128










\







.10 05 00 05 10
P 2



ao, 2.25; 22= 0.25; i= 128





--.. .-.. .


-10 -05 00 05 10


-10 05 00 05 LO
P12


,, =1.0;o22= 0.25; Kd= 1152


10 I05 00 05 L

Pl2


o,,= 2.25;o22 =0.25; = 1152



-


-10 -05 00 05 10

P,2


Legend: RESV ............ PB-20 PB-24 OA-24


Figure 5.4. Case 1 (C2) The Effect of Correlation on the Weighted Average MSE


S,,= 0,= 0.25; i = 128


-" -* - .






a22 = .25, and -.1 < P12 < .7 (top plot in left column), note that the OA24 performs

best. Also, the PB24 is in third place for all situations except when r,2 = 1152 and

all = U22 = .25 and all = 1.25, a22 = .25 (top and middle plots in right column) in

which case the PB20 outperforms it for large positive correlation. In Figure 5.4 where

the average MSE is multiplied by the number of design runs, the RESV16 performs

best in all situations with the OA24 performing second best. The PB20 performs

worst except for large positive correlation in which case the PB24 performs worst.

The superiority of the performance of the RESV16 in case 1 can be explained by

considering the alias matrix A*. For the RESV16, A* = 0. Consequently, it follows

from Section 5.1.1 that B1 = Q1 = 0, B2 = Q2 = 0, W2 = 0, and W3 = 0, where

the terms B1 and B2 are explained in Appendix B and W2 and W3 are explained

in Appendix C. As a result, the bias and mean squared error of Ty(xl) are equivalent

to the bias and mean squared error of (xa) obtained for the case in which the fitted

model was presumed to be correct. In other words, the second order control effects

present in the true model have no additional impact on the bias and mean squared

error of the estimated process variance when the RESV16 is utilized. For the OA24,

we also have A* = 0, and so the average bias and average MSE are equivalent to

that found in Chapter 3 for the correct model. However, the A* matrices for the

PB20 and PB24 designs are nonzero. This serves to explain the increase in bias and

MSE observed for these designs as compared to the case of the fitted model being

correct. It also explains why the OA24 performs better than the PB20 and PB24.

This differs from the results obtained in Chapter 3 where the OA24 was observed to

perform worst.

Figures 5.5 and 5.6 display the average bias for case 2 in which the true model

contains additional pure quadratic and interaction terms in the noise variables only.

In Figure 5.5, where the average bias is unweighted by the number of design runs, note

that the average biases are negative and also that a very small difference separates











S,,= 22= 0.25; 2= 1152


-10 -0- 5 00 0,5 10
S12


oi =1.0;22 = 0.25; K= 128


-1.0 -05 00 05 10
P 1.2


oa = 2.25; 0a22= 0.25; K 2= 128


-05 00

P 12


05 10


-L0 .0,5 00 05 1.0
P 12


oa = 1.0; 22= 0.25; K = 1152


-10 05O 0,0

P 12


01= 2.25;22 = 0.25; K = 1152


-10 -05 00

P 12


05 10


0.5 L


Legend: RESV ............ PB-20 ------ PB-24 OA-24


Figure 5.5. Case 2 (N2) The Effect of Correlation on the Average Bias


o = 022= 0.25; 2= 128









oa= 022=0.25; Ki= 128





.- -- -- -- -





/ \

-10 -0.5 0.0 05 1.0
P 12

11 =1.0;22 = 0.25; C2= 128


-10 0,5 00 05 1.0
S12


1 = 2.25; 22= 0.25; K= 128


0 -05 00 0.5 10


Legend: RESV ............


cH = 02= 0.25; K = 1152



..--- -- -- -- -- -- -









-10 05 00 0.5 1.0
P12

1 = 1.0; o = 0.25; K= 1152
/ ".




















10 -05 0.0 0.5 10
Pl2


S=2.25; a 22= 0.25; K2= 1152
a, = 2.25; a22 = 0.25; K< 2 = 1152


10 -.05 0.0 05 10


PB-20 -.-----. PB-24 OA-24


Figure 5.6. Case 2 (N2) The Effect of Correlation on the Weighted Average Bias


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







the best and worst designs. However, after multiplying by the number of design runs,

one can observe in Figure 5.6 that the RESV16 performs best followed by the PB20.

The OA24 and PB24 have nearly equal weighted average biases.

Figures 5.7 and 5.8 display the average MSE for case 2. In Figure 5.7 where the

average MSE is unweighted by the number of design runs, note that a very small

difference separates the best and worst designs although it appears that the Plackett-

Burman designs are performing better then the RESV16 and OA24. However, after

multiplying by the number of design runs, Figure 5.8 shows that the RESV16 has the

smallest average MSEs followed by the PB20 and PB24 and the OA24 comes in last.

In case 2 we observed again that A* = 0 for the RESV16 and OA24. Consequently,

it follows from Section 5.2.2 that Bi = Q, = 0 and B2 = Q2 Q5 = -Q5, where Q5

is positive definite and not dependent upon the design. Hence the additional impact

of the presence of second order noise effects in the true model is influenced by Q5.

The A* matrices for the PB20 and PB24 are nonzero. However, A* contains only

one nonzero column and so most of the elements of the matrices Q1 and Q2 are zero,

which serves to explain the negative biases that were observed for all the designs

considered.

Figures 5.9 and 5.10 display the average bias for case 3 in which the true model

contains additional interactions between noise factors and second order control effects.

In Figure 5.9, where the average bias is unweighted by the number of design runs,

one can observe that the PB24 performs best except for large negative correlation in

which case the PB20 performs best. The OA24 performs uniformly worst. However,

after multiplying by the number of design runs, Figure 5.10 shows that the RESV16

performs uniformly best and the OA24 performs uniformly worst.

Figures 5.11 and 5.12 display the average MSE for case 3. In Figure 5.11, where

the average MSE is unweighted by the number of design runs, we see that again

the PB24 performs best except for large negative correlation in which case the PB20










o =o 22,= 0.25; 2= 1152


-10 -05 0,0 0.5 10
P2 12


11 =1.0; o22 = 0.25; K 2= 128


12






















12

~


-10 -0,5 00 0.5 10
P t2


a = 2.25;a022= 0.25; K 2= 128


-1.0 -05 00

P 2


05 1.0


Legend: RESV


1.0 -0.5 00 05 10
P 12


o = 1.0;a22= 0.25; K 2= 1152


-1.0 -0.5 00 05 1

P 12


S= 2.25; 22 = 0.25; 2= 1152


-LO -05 0,0 05 10


............ PB-20 ------ PB-24 OA-24


Figure 5.7. Case 2 (N2) The Effect of Correlation on the Average MSE


(o = 2,= 0.25; K 2= 128










a = o22= 0.25; K 128
















P 12


0 = 1.0; a = 0.25; K 2= 128
\ /'
'\ 5














P t2


a = 2.25; 22= 0.25; K= 128

\- 0













-10 .05 00 05 10


0 =022= 0.25; K = 1152















-Lo 45 00 05 10
P12


1.0 = 0.25; K 1152
'/




























.10 05 00 05 10
P12


oG = 12.25; 22= 0.25; K= 1152















-10 405 00 05 10


Legend: RESV ............ PB-20 --- PB-2 OA-24


Figure 5.8. Case 2 (N2) The Effect of Correlation on the Weighted Average MSE











1,,= Oa= 0.25; i2= 128


















-10 -05 00 05 10
p12


a- = 1.0;a- =0.25;K2= 128



















-1.0 -0,5 0.0 0.5 1,0
P 12



ao, =2.250;o= 0.25; K 2= 128







-LO -05 0-0 0,5 1.0












10 .0.5 00 05 1)


o =o=22 0.25; K2= 1152










.- .
--










-10 -05 0.0 0.5 1.0

P 12


a = 1.0; a 22= 0.25; = 1152


-10 -05 0.0 0,5 1.0

P 12


a =2.25; 22 = 0.25; i = 1152


















-10 05 0.0 05 10
.10 .1)5 0.0 0.5 10


Legend: RESV ............ PB-20 ------- PB-24 OA-24


Figure 5.9. Case 3 (C2 x N) The Effect of Correlation on the Average Bias











1 =( =22 0.25; C 2= 128


















10 -05 0-0 05 10O
P 12


0, =1.0;22 = 0.25; K2= 128


-10 -05 00 0,5 10
P 2



a, = 2.25; a22= 0.25; K 2= 128


a1 =a22= 0.25; K = 1152


10 405 00 05 10
S12


C,, = 1.0; 022= 0.25; K = 1152


10 4)5 00 05 10

p,2


0, = 2.25; a22 = 0.25; K 2= 1152


00 0-5 10


-1.0 -0.5 0.0 05 10


r 12 r 12

Legend: RESV ............ PB-20 ----- PB-24 _- OA-24


Figure 5.10. Case 3 (C2 x N) The Effect of Correlation on the Weighted Average
Bias


-1.0 -05


-----------------------------

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












11, =a o= 0.25; K 2= 128
\
\
\













-1,0 -05 0.0 0,5 1,0
P 12


011 =1.0; a = 0.25; K2= 128

\
















1-0 05 00 0.5 1.0




a I= 2.25; a22= 0.25; I 2= 128
\

























10 05 00 05 10
N\
N-~


10 -0.5 00 0.5


10


a = 22= 0.25; K2= 1152


-10 0.5 00 05 1,0
P 12


Ca1 = 1.0;a = 0.25; K 2= 1152

\
\















-1,0 -0.5 00 0.5 1.0




(F= 2.25; 022 = 0.25; K 2= 1152


















-1'0 -0.5 0.0 0.5 10


Legend: RESV ............ PB-20 PB-24 OA-24


Figure 5.11. Case 3 (C2 x N) The Effect of Correlation on the Average MSE


N
N
N










o, = oa= 0.25; K2= 128


-10 -05 00 015 10
P12


0,, = 1.0; =0.25; K 2= 128


10 -05 00 05 10
P 12


11 = 2.25; a22= 0.25; K2= 128


-1.0 405 0.0 0-5 10


Legend: RESV ............


,=1 22= 0.25; K= 1152

\
\

\


10 05 00 05 1.0
S12


0,, =1.0; a = 0.25; K2= 1152

\

\.


.10 -05 00 0(5 1 I
P12


01 = 2.25; 022 = 0.25; K2= 1152


-10 -05 00 05 10


PB-20 ------ PB-24 OA-24


Figure 5.12. Case 3 (C2 x N) The Effect of Correlation on the Weighted Average
MSE


N

'N.
'N.







performs best. The OA24 performs uniformly worst. However, after multiplying by

the number of design runs, Figure 5.12 shows that the RESV16 performs best except

for large negative correlation in which case the PB20 performs best. The OA24

performs uniformly worst.

In case 3, all of the A* matrices were nonzero. The overall conclusions regarding

design performance are close to the recommendations obtained for the correct model

case with the exception of observing that the PB20 performs best for large negative

correlation. This was the first time we observed the PB20 to be the best design.

Figures 5.13 and 5.14 display the average bias for case 4 in which the true model

contains additional interactions between control factors and second order noise ef-

fects. In Figure 5.13, where the average bias is unweighted by the number of design

runs, one can observe that the average biases are negative and the PB20 performs

best except for large negative correlation in which case the PB24 performs best. The

RESV16 and OA24 average biases are nearly equal and these designs perform uni-

formly worst. However, after multiplying by the number of design runs, Figure 5.14

shows that the RESV16 performs best with two exceptions. The top plots in the left

and right columns show that the PB20 performs better than the RESV16 when the

noise variables are positively correlated. The OA24 performs uniformly worst.

Figures 5.15 and 5.16 display the average MSE for case 4. In Figure 5.15, where

the average MSE is unweighted by the number of design runs, observe again the PB20

performs best except for large negative correlation in which case the PB24 performs

best. The RESV16 and OA24 average MSEs are nearly equal and these designs

perform uniformly worst. However, after multiplying by the number of design runs,

Figure 5.16 shows that the RESV16 performs best with two exceptions. The top

plots in the left and right columns where oll = .25 show that the PB20 performs

better than the RESV16 except for large negative correlation. The OA24 performs

uniformly worst.