• TABLE OF CONTENTS
HIDE
 Title Page
 Acknowledgement
 Table of Contents
 Abstract
 Introduction
 Literature review -- testing for...
 An optimal check point method for...
 Use of near neighbor observations...
 Conclusions and recommendation...
 Appendices
 References
 Biographical sketch














Title: Testing lack of fit in a mixture model /
CITATION PDF VIEWER THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00097439/00001
 Material Information
Title: Testing lack of fit in a mixture model /
Physical Description: viii, 202 leaves : ill. ; 28 cm.
Language: English
Creator: Shelton, John Thomas, 1952-
Publication Date: 1982
Copyright Date: 1982
 Subjects
Subject: Mixtures   ( lcsh )
Mixtures -- Mathematical models   ( lcsh )
Statistics thesis Ph. D
Dissertations, Academic -- Statistics -- UF
Genre: bibliography   ( marcgt )
non-fiction   ( marcgt )
 Notes
Thesis: Thesis (Ph. D.)--University of Florida, 1982.
Bibliography: Bibliography: leaves 198-201.
Additional Physical Form: Also available on World Wide Web
Statement of Responsibility: by John Thomas Shelton.
General Note: Typescript.
General Note: Vita.
 Record Information
Bibliographic ID: UF00097439
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: alephbibnum - 000318861
oclc - 09205562
notis - ABU5710

Downloads

This item has the following downloads:

testinglackoffit00shel ( PDF )


Table of Contents
    Title Page
        Page i
        Page ii
    Acknowledgement
        Page iii
    Table of Contents
        Page iv
        Page v
        Page vi
    Abstract
        Page vii
        Page viii
    Introduction
        Page 1
        Page 2
        Page 3
        Page 4
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
        Page 18
    Literature review -- testing for lack of fit
        Page 19
        Page 20
        Page 21
        Page 22
        Page 23
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
        Page 29
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
        Page 36
        Page 37
        Page 38
        Page 39
    An optimal check point method for testing lack of fit in a mixture model
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
        Page 48
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
        Page 55
        Page 56
        Page 57
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
        Page 64
        Page 65
        Page 66
        Page 67
        Page 68
        Page 69
        Page 70
        Page 71
        Page 72
        Page 73
        Page 74
        Page 75
        Page 76
        Page 77
        Page 78
        Page 79
        Page 80
        Page 81
        Page 82
        Page 83
        Page 84
        Page 85
        Page 86
        Page 87
        Page 88
        Page 89
        Page 90
        Page 91
        Page 92
        Page 93
        Page 94
        Page 95
        Page 96
        Page 97
        Page 98
    Use of near neighbor observations for testing lack of fit
        Page 99
        Page 100
        Page 101
        Page 102
        Page 103
        Page 104
        Page 105
        Page 106
        Page 107
        Page 108
        Page 109
        Page 110
        Page 111
        Page 112
        Page 113
        Page 114
        Page 115
        Page 116
        Page 117
        Page 118
        Page 119
        Page 120
        Page 121
        Page 122
        Page 123
        Page 124
        Page 125
        Page 126
        Page 127
        Page 128
        Page 129
        Page 130
        Page 131
        Page 132
        Page 133
        Page 134
        Page 135
        Page 136
        Page 137
        Page 138
        Page 139
        Page 140
        Page 141
        Page 142
        Page 143
        Page 144
    Conclusions and recommendations
        Page 145
        Page 146
        Page 147
        Page 148
        Page 149
        Page 150
        Page 151
        Page 152
        Page 153
        Page 154
        Page 155
    Appendices
        Page 156
        Page 157
        Page 158
        Page 159
        Page 160
        Page 161
        Page 162
        Page 163
        Page 164
        Page 165
        Page 166
        Page 167
        Page 168
        Page 169
        Page 170
        Page 171
        Page 172
        Page 173
        Page 174
        Page 175
        Page 176
        Page 177
        Page 178
        Page 179
        Page 180
        Page 181
        Page 182
        Page 183
        Page 184
        Page 185
        Page 186
        Page 187
        Page 188
        Page 189
        Page 190
        Page 191
        Page 192
        Page 193
        Page 194
        Page 195
        Page 196
        Page 197
    References
        Page 198
        Page 199
        Page 200
        Page 201
    Biographical sketch
        Page 202
        Page 203
        Page 204
Full Text












TESTING LACK OF FIT IN A MIXTURE MODEL


BY

JOHN THOMAS SHELTON



















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


UNIVERSITY OF FLORIDA


1982































To Nydra

and

My Parents













ACKNOWLEDGEMENTS

I would like to express my deepest appreciation to Drs.

Andre' Khuri and John Cornell for suggesting this topic to

me and for providing constant guidance and assistance. They

have made this project not only a rewarding educational

experience but an enjoyable one as well. A special word of

thanks goes to Mrs. Carol Rozear for her diligence in

transforming my handwritten draft into an expertly typed

manuscript.


iii













TABLE OF CONTENTS
page

ACKNOWLEDGEMENTS............ ............. .. .............iii

ABSTRACT .............................................. vii

CHAPTER

ONE INTRODUCTION................................... 1

1.1 The Response Surface Problem............... 1
1.2 The Mixture Problem....................... 5
1.2.1 Mixture Models..................... 6
1.2.2 Experimental Designs for Mixtures.. 12
1.3 The Purpose of this Research:
Investigation of Procedures for Testing
a Model Fitted in a Mixture System for
Lack of Fit ............................... 17

TWO LITERATURE REVIEW--TESTING FOR LACK OF FIT..... 19

2.1 Introduction............................... 19
2.2 Partitioning the Residual Sum of Squares.. 21
2.3 Testing for Lack of Fit Without
Replicated Observations--Near Neighbor
Procedures ................................ 26
2.4 Testing for Lack of Fit with Check Points. 33

THREE AN OPTIMAL CHECK POINT METHOD FOR TESTING
LACK OF FIT IN A MIXTURE MODEL................. 40

3.1 Introduction............................... 40
3.2 Testing for Lack of Fit in the Presence
of an External Estimate of Experimental
Error Variation............................ 41
3.2.1 The Test Statistic................. 41
3.2.2 The Testing Procedure and an
Expression for the Power of
The Test........................... 45
3.2.3 A Method for Locating Optimal
Check Points ....................... 47
3.3 Testing for Lack of Fit When MSE Is
Used to Estimate Experimental Error
Variation..................................... 51
3.3.1 The Test Statistic................. 51
3.3.2 The Rejection Region and its
Relation to the Power of the Test.. 53
iv








3.3.3 A Method for Locating Optimal
Check Points ....................... 56
3.3.4 Determining Whether the Test Is
Upper Tailed or Lower Tailed....... 58
3.4 Examples. ................................ 67
3.4.1 Theoretical Examples............... 67
3.4.2 Numerical Examples................. 83
3.5 Discussion. ............................... 95

FOUR USE OF NEAR NEIGHBOR OBSERVATIONS FOR
TESTING LACK OF FIT .............................. 99

4.1 Introduction ............................. 99
4.2 Notation......... ..........................101
4.3 Shillington's Procedure ...................106
4.3.1 Development of MSEB ................109
4.3.2 Development of MSEW................110
4.4 Development of SSE(weighted) .............112
4.5 Equivalence of SSEW and SSEW(weighted)....116
4.6 The Test Statistic .........................118
4.7 The Testing Procedure and its Power........122
4.8 Selection of Near Neighbor Groupings....... 125
4.8.1 Example 1--Stack Loss Data.........129
4.8.2 Example 2--Glass Leaching Data.....134
4.9 Discussion................................ 142

FIVE CONCLUSIONS AND RECOMMENDATIONS ................145

APPENDICES

1 INFLUENCE OF X1 ON
P{F" > F } ...................156
v 1'2 ;11 '12 a;v ,v2

2 A CONTROLLED RANDOM SEARCH PROCEDURE FOR
GLOBAL OPTIMIZATION ............................159

3 STATISTICAL INDEPENDENCE OF d'V d/o2
2 -0-
AND SSE/o .................................... 164

4 THEOREM 3.1................................... 168

5 THEOREM 3.2...................................... 169

6 AN APPROXIMATION TO THE DOUBLY NONCENTRAL
F DISTRIBUTION.................................. 171

7 EQUIVALENCE OF SSEB AND SSLOF WHEN
REPLICATES REPLACE NEAR NEIGHBOR
OBSERVATIONS...................................... 172








8 LEMMA 4.1...................................... 175

9 PROOF OF THEOREM 4.1(i) ............ ......... 178

10 PROOF OF THEOREM 4.1(ii)......................182

11 PROOF OF THEOREM 4.1(iii).............. ........185

12 PROOF OF THEOREM 4.2.............................191
-1
13 PROOF OF THE EQUALITY SSE = d'V d + SSE .....193
A -0 -
REFERENCES .... ....................................... 198

BIOGRAPHICAL SKETCH.. .......... .. ..... ............202













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



TESTING LACK OF FIT IN A MIXTURE MODEL

By

John Thomas Shelton

May 1982



Chairman: Andre' I. Khuri
Cochairman: John A. Cornell
Major Department: Statistics

A common problem in modeling the response surface in

most systems, and in particular in a mixture system, is that

of detecting lack of fit, or inadequancy, of a fitted model

of the form E(Y) = X 1 in comparison to a model of the form

E(Y) = X1 + X 22 postulated as the true model. One method

for detecting lack of fit involves comparing the value of

the response observed at certain locations in the factor

space, called "check points," with the value of the response

that the fitted model predicts at these same check points.

The observations at the check points are used only for

testing lack of fit and are not used in fitting the model.

It is shown that under the usual assumptions of

independent and normally distributed errors, the lack of fit

test statistic which uses the data at the check points is an
vii







F statistic. When no lack of fit is present the statistic

possesses a central F distribution, but in general, in the

presence of lack of fit, the statistic possesses a doubly

noncentral F distribution. The power of this F test depends

on the location of the check points in the factor space

through its noncentrality parameters. A method of selecting

check points that maximize the power of the test for lack of

fit through their influence on the numerator noncentrality

parameter is developed.

A second method for detecting lack of fit relies on

replicated response observations. The residual sum of

squares from the fitted model is partitioned into a pure

error variation component and into a lack of fit variation

component. Lack of fit is detected if the lack of fit

variation is large in comparison to the pure error

variation. This method can be generalized when "near

neighbor" observations must be substituted for replicates.

In this case, the test statistic (assuming independent and

normally distributed errors) has a central F distribution

when the fitted model is adequate and a doubly noncentral F

distribution under lack of fit. The arrangement of near

neighbors is seen to affect the testing procedure and its

power.


viii














CHAPTER ONE
INTRODUCTION

1.1 The Response Surface Problem

A mixture problem is a special type of a response

surface problem. First we shall define the general response

surface problem and indicate the basic objectives sought in

its analysis, and follow this development with a discussion

of the mixture problem.

In the general response surface problem, we are inter-

ested in studying the relationship between an observable

response, Y, and a set of q independent variables or

factors, xl, x2, ..., Xq, whose levels are assumed con-

trolled by the experimenter. The independent variables are

quantitative and continuous. We express this relationship

in terms of a continuous response function, p, as



Y = (ul, x2 uq ) + u



where Yu is the uth of N observations of the response col-

lected in an experiment, and xui represents the uth level of

the ith independent variable, u = 1, 2, ..., N, i = 1, 2,

..., q. The exact functional relationship, p, is unknown.

The term Eu is the experimental error of the uth







observation. It is assumed that E(cu) = 0, E(Eueu,) = 0,

for u u', and E(e ) = 2, for u = 1, 2, ..., N.

As the form of ( is unknown and may be quite complex, a

low order polynomial (usually first or second order) in the

independent variables xl, x2, ..., Xq is generally used to

approximate p. This may be justified by noting that such

polynomials constitute low order terms of a Taylor series

expansion of < about the point xl = x2 = ... = xq = 0,

(Myers, 1971, p. 62). Cochran and Cox (1957, p. 336) point

out that these low order polynomials may give a poor approx-

imation to q when extrapolated beyond the experimental

region, and thus should not be used for this purpose.

A linear response surface model may be written in

matrix notation as



Y = X4 + E (1.1)



where Y is an Nxl vector of observable response values, X is

an Nxp matrix of known constants, a is a pxl vector of

unknown parameters (regression coefficients), and a is the

Nxl vector of random errors. When the model is a first or a

second degree polynomial, the columns of X correspond to the

first or second degree powers of the independent variables

xl, x2, ..., Xq, or their cross products. If the model

contains a constant term, 80, the first column of X will

correspond to this term, and will consist of N ones. Since

E(c) = 0, an alternative representation for the response







surface model of (1.1) is



E(Y) = X .



Once the form of the model that will be used to approx-

imate ((xl, x2, ..., Xq) is chosen, the next step is to

estimate the regression coefficients, a, and then use the

estimated model to make inferences about the true response

function, (. The estimation of the elements of a is usually

accomplished by ordinary least squares techniques. For the

purpose of testing hypotheses concerning the regression

coefficients, ., it is assumed that L has a normal distribu-

tion, that is, e ~ NN(O, a 2N).

Perhaps the most common objective in the exploration of

a response system is the determination of its optimum

operating conditions. By this we mean that it is desired to

find the settings of xl, x2, ..., xq that optimize (, which

in some applications may be interpreted as maximizing p,

while in other applications a minimum value of ( may be of

interest. It is also often desirable to determine the be-

havior of the response function in the neighborhood of the

optimum. For second order response models, such an investi-

gation can be carried out by performing a canonical analysis

of the second order surface as discussed in Myers (1971).

For simple systems having only one or two independent

variables, the response surface may be explored by just

plotting the fitted response values against values taken by







the independent variables. If q = 1, implying only one

independent variable, say xl, then a two-dimensional plot of

the fitted response against xl may be used to locate the

optimum, as well as to investigate the response behavior in

other parts of the experimental range of xl. If q = 2, and

the two independent variables are xl and x2, then a plot of

the contours of constant response over the region specified

by the ranges of the values for xI and x2 can be used to

describe the response surface.

The properties that the fitted model possesses in terms

of its ability to represent the true surface, p, depend on

the settings of xl, x2, ..., Xq at which values of Y are

observed. Thus the experimental design is of great impor-

tance. Much work has been done on the construction of

designs that are optimal with respect to one criterion or

another involving the fitted response and/or the true unfit-

ted model. Box and Draper (1975) list fourteen criteria to

consider when choosing a design for investigating response

surfaces. Myers (1971) gives several designs for fitting

first and second order polynomial models. A discussion of

specific design considerations will not be attempted here,

as such a discussion is not the focus of this dissertation,

and would necessarily be lengthy.

The initial steps in the analysis of a response system

may be described as follows: First an attempt is made to

approximate the true response function, p(xl, x2, ..., Xq),

usually with a low order polynomial in xl, x2, ..., xq.







After the form of the model has been chosen, then comes the

selection of an appropriate experimental design, which

specifies the settings of the independent variables at which

observed values of the response will be collected. The

observed values of the response are used in estimating the

regression coefficients in the model, using, in general,

ordinary least squares. After a test for "goodness of fit"

of the model verifies the fitted model is adequate, the

fitted model is used in determining optimum operating condi-

tions for the response system.

1.2 The Mixture Problem

A mixture system is a response system in which the

response depends only on the relative proportions of the

components or ingredients present in a mixture, and not on

the total amount of the mixture. For example, the response

might be the octane rating of a blend of gasolines where the

rating is a function only of the relative percentages of the

gasoline types present in the blend. The proportion of each

ingredient in the mixture, denoted by xi, must lie between

zero and unity, i = 1, 2, ..., q. The sum of the propor-

tions of all the components will equal unity, that is,


q
0 < x. < 1, i = 1,2,...,q, E x. = 1. (1.2)
i=l


The factor space containing the q components is represented

by a (q 1)-dimensional simplex. For q = 2 components, the

factor space is a straight line, whereas for q = 3







components, the factor space is an equilateral triangle, and

for q = 4 components, the factor space is represented by a

regular tetrahedron.

The objectives in the analysis of a mixture response

system are, in general, the same as in any response surface

exploration. That is, one seeks to approximate the surface

with a model equation by fitting an equation to observations

taken at preselected combinations of the mixture com-

ponents. Another objective is to determine the roles played

by the individual components. We shall not concern our-

selves with this but rather concentrate on the empirical

model fit. Once the model equation is deemed adequate an

attempt is made to determine which of the component combina-

tions yield the optimal response. The models used to repre-

sent the response in a mixture system are in most cases

different in form from the standard polynomial models. The

first type of model form that we discuss is the canonical

polynomial suggested by Scheffe.

1.2.1 Mixture Models

Scheffe (1958) introduced a canonical form of the poly-

nomial model for representing the response in a mixture

system. These canonical polynomial models are derived from

the standard polynomials using the restrictions on the xi

shown in (1.2). With q = 2 mixture components, for example,

the standard second order polynomial model is of the form


2 2
E(Y) = a0 + alx1 + a2x2 + a12 x2 + all x + a22x2 (1.3)







Restrictions (1.2) imply that a0 = a0(xl + x2)'

x2 = l(- and x2 = 2(- thus (1.3) can be
1 xl(l x2), x2 x2(l -xl)
written in the canonical form



E(Y) = 81x1 + 2x2 + 12XlX2



where 1= a0 + al + all' 82 = a0 + a2 + a22' and 812= al2

- a11 a22. There is no constant term in the above canoni-

cal form and the pure quadratic terms in equation (1.3) have

been absorbed in the xixj terms.

The general form of the canonical polynomial of degree

d in q mixture components can be written as


q
E(Y) = E 8ixi, for d = 1,
i=l

q q
E(Y) = Z i.x. + Z E .x.x. for d = 2, and
i=l 1. i
q q q
E(Y) = a .x. + EZ .x.x. + E E 6..x.x.(x. x.)
i=l 1 i
q
+ E a ijkx.x.xk for d = 3. (1.4)
1Ki

The fourth degree canonical polynomial in q components is

given in Cornell (1981, p. 64). The general canonical poly-

nomial of degree d > 4 in q components does not explicitly

appear in the literature, but is mentioned in Scheffe

(1958). If terms of the form 6ijxixj(xi xj) are removed

from the full cubic model (1.4), then the remaining terms







make up what is referred to as the special cubic model. For

example, for q = 3 components, the special cubic model is



E(Y) = 1x1 + 82x2 + 3x3 + a12 2 + 813X1X3



+ 823x2x3 + 123XlX2X3



Scheffe's canonical polynomial models are used for

approximating the response surface in many mixture systems.

Their popularity stems from the ease in interpreting the

coefficient estimates, especially when the models are fitted

to data collected at the points of the associated designs

(see Section 1.2.2). However, other models have been intro-

duced which seem to better represent the response when the

components have strictly additive blending effects. We

present some of them now.

Becker (1968) introduced three forms of homogeneous

models of degree one which he recommends instead of the

polynomial models when one or more of the mixture components

have an additive effect or when one or more components are

inert. A function f(x, y, ..., z) is said to be homogeneous

of degree n when f(tx, ty, ..., tz) = tnf(x, y, ..., z), for

every positive value of t and (x, y, ..., z) (0, 0, ...,

0). These models, which Becker refers to as models HI, H2,

and H3, are of the form








q q
HI: E(Y) = Z B.x. + Z Bijmin(xi, x ) + ...
i=1 1 i

+ 12...qmin(xl, x 2, ..., x )


q q
q q 2-1
H2: E(Y) = E Bx. + Z Bix. .x ./(x. + x .) +
i= 13 1 3 1
i=1 14i

+ 812...qxlx2..x /(x + x2 + ... + xq)q-


q q 1
H3: E(Y) = 6 .x. + Z . (x.x 1/2
i=1 l i

+ Bl2...q(X l2...xq)/q


Each term in the H2 model is defined to be zero when the

denominator of the term is zero.

Draper and St. John (1977) suggest a model which in-

cludes inverse terms, 1/xi, in addition to terms in the

Scheffe polynomials. Such a term is used to model an

extreme change in the response as xi approaches zero. The

experimental region of interest is assumed to include the

region near the zero boundary (xi = 0), but does not include

the boundary itself. One example of this type of model is

the Scheffe linear polynomial model with inverse terms


q q -1
E(Y) = Z Bixi + Z B .x .
i=l i=l -







Another model form that is useful in the study of the

response in a mixture system is the model containing ratios

of the component proportions. A term such as xi/xj measures

the relationship of xi to xj rather than the percentage of

each in the blends. Snee (1973) points out that the ratio

model presents a useful alternative to the Scheffe and

Becker models in that the ratio model describes a different

type of curvature. He notes that the curvilinear terms for

the Scheffe and Becker models, when plotted as a function of

xi, are symmetric functions about xi = 1/2, whereas the

ratio term xi/xj is a monotone function when plotted against

xi.

The terms in the ratio models may also contain sums of

the components. For example, with q = 3 components, we

might express the second order model


q-1 q-1 q-1
E(Y) = B + E izi + EE ijziz. + E 68iz.
i=l l

(note the constant term) where zl and z2 are defined as

z1 = xl/(x2 + x3) and z2 = x2/x3. Some terms will be unde-

fined if points from the boundary of the experimental sim-

plex are included in the design, for example, if x3 = 0,

then z2 = x2/x3 is not defined. Snee (1973) suggests adding

a small positive quantity, c, to each xi in this case.

This, of course, will not be of concern if the experimental

region is entirely inside of the simplex.







When one or more of the components is inactive, Becker

(1978) suggests that a ratio model that is homogeneous of

degree zero in the remaining components is appropriate. In

three components, such a model is of the form



E(Y) = 80 + 1X/(x1 + x2) + 2x2/(x2 + x3)


3
+ 83x /(x + x ) + Z BZ ..h..(x., x )
s3 3 1 3 1i li

+ 8123h123(xl, x2, x3), (1.5)



where hij and h123 are specified functions that are homoge-

neous of degree zero. The function h123 is intended to

represent the joint effect of all three components simulta-

neously. If in fitting a model of the form (1.5) we deter-

mine the model should be



E(Y) = 80 + a1X1/(X1 + x2) + 812h12(xl, x2)



then component three is said to be inactive and is removed

from further consideration. The model of equation (1.5) may

produce an extreme value near the vertices of the simplex

factor space when there are no inactive components. In this

case it is suggested that a model of the form (1.5) be used

only when the proportions are restricted so that no two of

the xi are simultaneously very close to zero. Becker notes

that other authors who have suggested ratio models have also







used them primarily over a subregion inside the simplex

factor space. Apparently this is where they are most appro-

priate.

1.2.2 Experimental Designs for Mixtures

As in the general response surface problem, one of the

major concerns in exploring a mixture system is that of

choosing the experimental design for collecting observed

values of the response that will be used in fitting the

model. Scheffe (1958) proposed the {q,m} simplex lattice

designs for exploring the entire q-component simplex factor

space. In these designs, the proportions used for each

component have the m + 1 values spaced equally from zero to

one, xi = 0, 1/m, 2/m, ..., (m l)/m, 1, and all possible

mixtures with these proportions for each component are

used. The number of design points in the {q,m} simplex

lattice design is (m + q 1). The main appeal of these
m
designs is that they provide a uniform coverage of the fac-

tor space. Another feature, which Scheffe (1958) demon-

strates, is that the parameters of the canonical polynomial

of degree m in q components are expressible as simple linear

combinations of the true response values at the design

points of the {q,m} simplex lattice. The {3,2} simplex

lattice, which consists of six design points, is represented

on the two dimensional simplex in Figure 1 along with the

triangular coordinates (xl, x2, x3).

Scheffe (1963) also developed the simplex centroid

designs consisting of 2q 1 points, where the only mixtures







considered are the ones in which the components present

appear in equal proportions. That is, in a q-component

simplex centroid design, the design points correspond to the
q
q permutations of (1, 0, 0, ..., 0), the () permutations of

(1/2, 1/2, 0, ..., 0), the (3) permutations of (1/3, 1/3,

1/3, 0, ..., 0), ..., and the point (1/q, l/q, ..., l/q).

This design alleviates the problem inherent in the {q,m}

simplex lattice designs of observing responses at mixtures

containing at most m components. To give an example, the

q = 3 simplex centroid design is made up of 23 1 = 7

design points, and is equivalent to the {3,2} simplex

lattice design augmented by the center point (xl, x2, x3)

(1/3, 1/3, 1/3). This design is represented in Figure 2.

Scheffe (1963) mentions that the number of parameters

in the polynomial model of the form


q q q
E(Y) = .ix. + EE B. .x.x. + EE 8ikx. xxk
i=l 1 i

+ ... + B12...qxlx2 ... q (1.6)



is 2q 1 and therefore these models have a special rela-

tionship with the simplex centroid design in q components.

This relationship is that the number of terms in the model

equals the number of points in the design and as a result

the parameters in model (1.6) are expressible as simple

functions of the responses at the 2q 1 points of the sim-

plex centroid design. Polynomial models of the form (1.6)



















(2 2


(0,1,0)




Figure 1.













(





(0,1 ,0)


(0,0,1 )


x =I ( x, I )
2 2'2 3


The {3,2}


simplex lattice design.


x =I


(0,0,1)


X3


The q = 3 simplex centroid design.


x I)
~~2 2'2


Figure 2.







therefore are natural models to fit using the simplex cen-

troid design.

Ratio models may be desirable when the interest in one

or more of the mixture components is in terms of their rela-

tionship to one another, rather than in terms of their per-

centages in blends. Kenworthy (1963) proposed factorial

arrangements for ratio variables. An example of the use of

ratios is the following three component system in which the

mixture components are constrained by the upper and lower

bounds:



.2 < x1 < .4, .2 < x2 < .4, .3 < x3 < .5. (1.7)



The ratio variables of interest are zl = x2/xl and

z2 = x2/x3, and it is desired to fit either a first or a

second order polynomial model in zI and z2. For such a

problem, we can define a 22 and a 32 factorial design that

can be used for fitting the first and second order poly-

nomial models, respectively, by taking as design points the

intersection of rays passing from two of the three vertices

of the two-dimensional simplex through the region of

interest defined by the constraints (1.7). Kenworthy's 22

factorial design is shown in Figure 3.

Becker (1978) uses rays extending from one or more

vertices of the simplex factor space to the opposite bound-

aries in developing "radial designs." These designs are

suggested for detecting the presence of an inactive











x :




0 Design Points









x :1 x :I
2 3


Figure 3. Kenworthy's 22 factorial design.



component or in another case a component which has an addi-

tive effect, when models containing ratio terms that are

homogeneous of degree zero are fitted.

McLean and Anderson (1966) suggest an algorithm for

locating the vertices of a restricted region of the simplex

factor space which is defined by the placing of upper and

lower bounds on the mixture component proportions. The

vertices of the factor space and convex combinations of the

vertices are the candidates for design points for fitting a

first or second degree polynomial model in the mixture com-

ponents. One drawback of the "extreme vertices" design is

that the design points are not uniformly distributed over

the factor space resulting in an imbalance in the variances

of Y(x), see Cornell (1973).







Another method that has been suggested for studying the

response over a sub-region of the simplex mixture space is

to transform the q mixture components into q 1 independent

variables. Transforming to an independent variable system

was first suggested by Claringbold (1955) and later proposed

by Draper and Lawrence (1965a, 1965b) and Thompson and Myers

(1968). Standard response surface polynomial models in the

transformed variables can be fitted to data values collected

on standard designs and a design criterion such as the aver-

age mean square error of the response can be employed when

distinguishing between designs. Thompson and Myers (1968)

suggest the use of rotatable designs (see also Cornell and

Good, 1970).

Designs other than rotatable designs, such as multiple

lattices and symmetric-simplex designs, to name a few, have

been suggested in the literature for fitting models to a

mixture system which may be appropriate depending on par-

ticular experimental situations. However, as the intent

here is not to give an exhaustive list but only a sampling

of available designs, we shall not discuss designs further

but instead state the purpose of this work.

1.3 The Purpose of this Research:
Investigation of Procedures tor Testing a Model
Fitted in A Mixture System for Lack of Fit

A common problem in modeling the response in a mixture

system is that of detecting lack of fit, or inadequacy, of a

fitted model of the form E(Y) = XI when the true model is

of the form E(Y) = X2 + X2 2. The statistical literature







suggests several procedures for testing lack of fit, which

will be described in Chapter Two. In general, the authors

of these procedures are not specific in stating hypotheses

to be tested and do not adequately discuss the power of

their procedures.

The major purpose of this research is to investigate

the power of two of the testing procedures appearing in the

literature in detecting the inadequacy of a fitted model

when the general form of the true model is specified. Our

findings for a "check points" lack of fit testing procedure

are presented in Chapter Three while Chapter Four contains

findings for a "near neighbor" lack of fit testing proce-

dure. For both procedures, explicit formulas for the power

of the test are given in terms of cumulative probabilities

of either the noncentral F or doubly noncentral F distribu-

tion, which are derived by assuming that the response obser-

vations are independent and normally distributed. Addition-

ally, we propose methods for maximizing the power of the

testing procedures. In the final chapter, we make some

concluding comments concerning both of these procedures.













CHAPTER TWO
LITERATURE REVIEW--TESTING FOR LACK OF FIT

2.1 Introduction

Let us return to the general response surface problem

and assume the true response is to be approximated by

fitting a model of the form



E(Y) = XI1 (2.1)



where Y is an Nxl vector of observable response values, X is

an Nxp matrix of known constants, and El is a pxl vector of

unknown regression coefficients. We wish to consider the

situation in which the true model contains terms in addition

to those in the fitted model. Then the true model has the

form



E(Y) = X1 + X2a2 (2.2)



where X2 is an Nxp2 matrix of known constants, and 12 is a

P2xl vector of unknown regression coefficients. We assume
that the vector Y has the normal distribution with

var(Y) = a 2N

It is desirable to determine the suitability of the

fitted model given by Eq. (2.1) when in reality the true

model is of the form given by Eq. (2.2). The process of
19







making this determination is referred to as testing for lack

of fit of the fitted model.

There are three general approaches to testing for lack

of fit. The first approach requires that there be replicate

observations of the response at one or more design points,

and involves partitioning the residual sum of squares from

the fitted model into a sum of squares due to lack of fit

and a sum of squares due to pure error. A large value for

the ratio of the mean square due to lack of fit to the mean

square due to pure error provides evidence for lack of fit.

If replicate observations are not available then the

above approach to testing for lack of fit cannot be used.

Green (1971), Daniel and Wood (1971), and Shillington (1979)

have proposed alternative methods that are applicable in

this case. Their approach is to group values of the

response which are observed at similar settings of the

independent variables and to call these grouped values

"pseudoreplicates" or "near neighbor observations." They

then treat these pseudoreplicates as they would treat true

replicates to form statistics for lack of fit testing,

although arriving at their respective statistics through

different approaches.

A third approach to testing for lack of fit involves

the use of "check points." In this method a model of the

form (2.1) is fitted to data at the design points and

additional observations are collected at other points in the

experimental region. The additional points other than the







design points are called check points, and the data at these

check points are not used in fitting the model. Lack of fit

is tested by comparing the values of the response observed

at the check points to the values of the response which the

fitted model predicts at these same check points.

We now discuss the first method mentioned above of

testing for lack of fit which involves partitioning the

residual sum of squares.

2.2 Partitioning the Residual Sum of Squares

The method for testing lack of fit which makes use of a

partitioning of the residual sum of squares from the fitted

model requires there be replicate observations of the

response at some of the design points (Draper and Smith,

1981, p. 120). When a model of the form (2.1) is fitted,

the residual sum of squares is defined as

n.
n 1 2
SSE = EZ (Y Y )
i=1 j=l 1

-1
=Y'(I' X(X'X) X')Y



where n is the number of distinct design points, ni > 1 is

the number of replicate observations at the ith design

point, Yij is the jth observed value of the response at the

ith design point, Yi is the value which the model of the

form in Eq. (2.1), fitted by ordinary least squares

techniques, predicts for the response at the ith design
n
point, and N = Z n. Using the replicated observations
i=l 1







only, a pure error sum of squares can be calculated as


n "i 2
SSE = E E (Y. Y. )
pure i=l 13 '
i=1 j=1


where Yi. is the average of the values of the response

observed at the ith design point. The sum of squares due to

lack of fit can be obtained by taking the difference


LOF


=SSE SSE
pure


This partitioning of the residual sum of squares is

displayed in the analysis of variance table in Table 1.


Table 1. Analysis of Variance--
Partitioning the Residual Sum of Squares.


Source
of Variation

Regression

Residual

Pure Error

Lack of Fit

Total(corrected)


Sum
of Squares

b{X'Y (1'Y)2/N

SSE

SSEpure

SSLOF
Y'Y (l'Y)2/N


Degrees
of Freedom

p 1

N p

N n

n p

N- 1


bl represents the ordinary least squares estimator of B in
-1
model (2.1), b = (X'X) X'Y, and 1 is an Nxl vector of
ones.
ones.


Mean
Square



MSE

MSEpure

MSLOF







To test the hypothesis of zero lack of fit, that is

HO: lack of fit = 0 or E(X) = XEI, an F statistic is formed


MS
LOF
F = F (2.3)
MSE
pure


which possesses a central F distribution if the true model

is of the form (2.1), but has a noncentral F distribution if

the true model is of the form (2.2). In other words



F ~ F
n-p,N-n


under H : E(Y) = X8_ and



F ~ F'
n-p,N-n;X2


under H : E(Y) = XB + X2_ where X2 is the noncentrality
a 1+ 22 2
parameter 2 = B(X2-XA)'(X2-XA)B2/22 and A = (X'X)- X'X2

Under H E(MSLOF) = 2 + (X2 XA)'(X2 XA) 2/(n-p) and

E(MSEpure) = o2 (Draper and Smith, 1981, p. 120), hence HO

is rejected in favor of Ha if the value of F in (2.3)

exceeds the upper 100a percentage point of the central F

distribution, Fa;n-p,N-n. When HO is rejected, we conclude

that a significant lack of fit is present.

Draper and Herzberg (1971) demonstrated that the lack

of fit sum of squares can be partitioned into two

statistically independent sums of squares, SSL1 and SSL2'

when there are replicate observations at the center of the







response surface design and when the basic design without

center points is nonsingular. If the true model and the

fitted model are of the same form as in equation (2.1) then

the two F ratios FL1 = [SSLl/(n p 1)]/MSEpure and

FL2 = SSL2/MSEpure are both distributed as central F random

variates, with respective numerator and denominator degrees

of freedom (n p 1), (N n) for FL1 and 1, (N n) for

FL2. If the true model is of the form shown in equation

(2.2), then FL1 and FL2 are both distributed as noncentral F

random variates. The expected values of SSL1 and SSL2 are

used to show what functions of E2 are testable with FL1 and

FL2*

Two examples are presented by Draper and Herzberg to

illustrate this testing for lack of fit. The first example

makes use of a first order orthogonal design in k factors

augmented with center point replicates for fitting a first

order polynomial model. If the true model is of the second

order, then FL2 can be used to test a hypothesis concerning

the parameters associated with the pure quadratic terms in

the model. If all such parameters are zero, then FL1

provides a check on interaction terms. The second example

illustrates the fitting of a second order polynomial model

to a second order design with all odd design moments of

order six or less zero. If the true model is third degree,

then FL1 can be used to test the significance of the third

order terms, while FL2 tests terms of order greater than

three. The partitioning of SSLOF into SSL1 and SSL2 and the







corresponding tests of hypotheses are also given in Myers

(1971, p. 114-119), for the special case of fitting a first

order polynomial model to a 2q factorial or a fraction of a

2q factorial design augmented with center point replicates

and the true model is of the second degree.

A more complete partitioning of the lack of fit sum of

squares in an attempt to obtain a more detailed diagnosis of

the lack of fit of the fitted model is given in a technical

report written by Khuri and Cornell (1981). The lack of fit

sum of squares, which has n p degrees of freedom, is

partitioned into n p independent sums of squares, each

having one degree of freedom. The expected values of these

single degree-of-freedom sums of squares are used to

identify at most n p linearly independent causes for the

lack of fit variation. Tests of significance are performed

on the assumed contributing causes. This method enables the

screening of all subsets of E2 in order to identify those

subsets which are most responsible for lack of fit of the

fitted model.

We shall now discuss the second general approach used

in lack of fit testing, which is to test for lack of fit by

making use of response values observed at points which are

near neighbors in the factor space when true replicate

observations are not available.







2.3 Testing for Lack of Fit Without
Replicated Observations--Near Neighbor Procedures

Green (1971) suggests the following approach when

testing for lack of fit if there are no design points at

which replicate observations of the response are

available. The N observed values of a response, Y,

considered a function of only one variable, x, are divided

into g groups, by grouping observations which have similar

values of x. Green hypothesizes a model of the form Y= Ha +

e, where Y is an Nxl vector of observable responses, H is an

Nxm matrix whose columns correspond to known functions of

the variable, x, is an mxl vector of unknown regression

coefficients, and e is the Nxl vector of random errors,
2
N ~ N (0, o N).

Green's method assumes that the vector of differences

(EY Hg) can be well approximated by a dth order polynomial

in x within each of the g groups, d > 1. An alternative

model of the form



Y = H v + n +



is given, where is distributed as NN(Q, o21N), H1 is an

Nx [g(d + 1) + ml] matrix of known constants, u is a

[g(d + 1) + ml]xl vector of regression coefficients, and .,

as Green states is "a small vector." The first g(d + 1)

columns of H1 correspond to the polynomial terms for the g

groups (with (d + 1) terms for each group), the rightmost

mI < m columns in H1 correspond to terms that are in the







fitted model, but are not represented among the g(d + 1)

polynomial terms in the alternative model.

Under the assumption that a = Q, the presence of lack

of fit is tested by using the test statistic:




Y'[H (HHl) -H H(H'H) 1H']Y/[g(d + 1) + mi m]
F = *----------------------
-1
Y'[I H (HH1I) H{]Y/[N g(d + 1) mi]


(2.4)



This statistic is of the same form as the F statistic used

in the standard multiple regression test of a postulated

model against a more general one which includes the

postulated model as a special case. Lack of fit is

suspected if the calculated F ratio in (2.4) is greater

than Fa;g(d+l)+ml-m, N-g(d+l)-ml where this latter quantity

is the upper 100a percentage point of the central F

distribution.

Green notes that when there is no lack of fit, the

quadratic forms Y'[H1(HIH ) -H H(H'H) H']Y and

Y'[I H (HIH,) Hi]Y are distributed independently as

o2X2 with g(d + 1) + mi m and N g(d + 1) mi degrees of

freedom, respectively. In this case the F ratio in (2.4)

possesses a central F distribution. If there is lack of fit

on the other hand, then these two quadratic forms are

distributed as noncentral chi-squares, multiplied by 02,

with respective noncentrality parameters







I = [H + '[H (HiH1) -H' H(H'H)- H'][H + n]
and ; = n'[I H1(HH ) -HI]n Thus the assumption that

n = 0 can affect the power of the test, since if n # 0 the

expected value of MSE is greater than a2, where MSE is the

quadratic form in the denominator of the F ratio. Hence if

n 0 the probability of calculating a large F value is

reduced, and we are less likely to detect lack of fit using

an upper tailed rejection region.

Daniel and Wood (1971) suggest another method for lack

of fit testing when replicated observations of the response

are not available. They make use of "near replicates" to

obtain an estimate of a, which is the standard deviation of

the observable responses in the true model. The value of

the estimate a is compared to the square root of the

residual mean square from the analysis of the fitted

model. Lack of fit is indicated if the square root of the

residual mean square is large compared to the estimate o.

To determine when observations are near replicates so that

an estimate of a can be found, they define the squared

distance between any two data points, j and j', to be

measured by

K

D2 = x. ij.)/5y]2
3J' i=1


where xij and xij, are the values of the ith independent

variable corresponding to the observations yj and yj,,

respectively, i = 1, 2, ..., K, and bi is the ordinary least







squares estimate of the ith regression coefficient. In the

denominator, s is the square root of the residual mean

square for the fitted model.

To obtain an estimate of a from near replicates, let

And = |dj d ,j, n = 1, 2, ..., (2), where dj and d are

the residuals at points j and j', respectively, and where

there are N data observations in the experiment. Since the

expected value of the range for pairs of independent

observations from a normal distribution is 1.1280, a running

average of the And's is calculated and their average is

multiplied by .886 = (1/1.128) to get a running estimate,

sn, of a. That is, s = .886 n And/n The closest pair

of observations as judged by Dj, is used to begin the

running estimate, the next closest pair (next "nearest

neighbors") is used for A2d, and the procedure continues

until sn "stabilizes." The stabilized value of sn is used

to estimate a.

A third method for testing for lack of fit without

replication is given by Shillington (1979). The fitted

model is of the form

Y = XB + E (2.5)

where Y (Nxl), X (Nxp), and (pxl) are defined as in (1.2)

and E ~ N (0, a IN) The test for lack of fit of the

fitted model is a test for whether the true model has the

form


y = X1 + 6 + E 1







where 6 (Nxl) is a fixed effect quantifying the departure of

(2.5) from the true model.

Shillington assumes that the data can be grouped into g

cells, with nj observations in the jth cell, determined in

advance. Letting Cj refer to the jth cell, j = 1, 2, ...,

g, a vector of cell averages is written YC (gxl), where the

jth element of YC is the average of the observed responses

in Cj. The matrix XC of independent variables associated

with Y is the gxp matrix where the elements in the jth row
n.
are x'. = i x!./n. that is, row j of X is the row
i-. =1 -
vector x'. The matrix XC is assumed to be of full rank

p < g. Also within each cell are defined the differences

W.. = Y.. Y. i C. j = 1, 2, ..., g, where Y is

the jth element of Yc.

The two independent data sets, YC and {Wij} with g and

N g degrees of freedom, respectively, are used to find two

independent estimates of o2. The first estimate is written

as


9 2
MSEB = E n (Y x' /(g -)
j=l B


where g4 is the weighted least squares estimate of g using

the regression of cell means, YC, on XC. The second

estimate of o2 uses the within cell deviations on cell

means, {Wij}, and is

n.
g n 2
MSEJ = Z (W. W. ) /(N g r),
j=1 i=l







where r is the rank of an Nxp matrix with rows equal to

x! -x' i C C., j = 1, 2, .., g. If the matrix of
-13j -.j 3
independent variables, corrected for cell means, is of full

rank, then r = p. Here 7ij is the estimate of Wij from the

regression of cell residuals {Wij} on the associated vectors

of independent variates, x'. x'.

If the fitted model is the correct model, then MSEB and

MSEW are independent estimates of a2 and the ratio MSEB/MSEW

is an F statistic with g p and N g r degrees of

freedom. When all observations in a cell have the same

settings of the independent variables, that is, the

observations are truly replicates for all cells, then this F

statistic is identical to the F statistic in the usual lack

of fit test in which the residual sum of squares is

partitioned into lack of fit and pure error sums of squares,

as given in Draper and Smith (1981, p. 120).

If the true model is Y = X8 + 6 + E however, and if

we let X'6 = 0 and 2 = 6/N, then



E(MSEB) = o2 + [I XC(X'XC-)1'] B/(g-P)

n.
where B (gxl) has jth component equal to E 6../n .
-B i=
Furthermore, with this latter true model form



E(MSE ) = 02 + 67(I XW(XJ XW)- X)6W/(N g r)







where _W has the components ij .j, i e Cj, j = 1, 2,

..., g. The matrix XW (Nxp) has the rows x x'

i e C., j = 1, 2, ..., g. The power of the F test,

F = MSEB/MSEW, depends on the relative bias of the estimates

of 2, that is, the biases in MSEB and MSEW.

Shillington states that the power of the F test which

makes use of F = MSEB/MSEW is maximized by forming cells so

that the bias of E(MSEW) is minimized. This is the same as

forming cells so that the within cell variation in 6 is

minimized. Shillington (1979, p. 141) also states,

"Observations with near covariate (independent variable)

values might be expected to have similar 6 values, since we

assume that 6 varies in some continuous but unknown fashion

with X. This justifies the usual procedure of forming

groups by collapsing observations with adjacent covariate

values. Indeed, if covariates do not vary within cells we

have the usual lack of fit test and maximum power."

By imposing a further structure on the form of 6, it is

shown that if the F test has an upper tailed rejection

region, the power is maximized by selecting the group sizes

as nj = 2, j = 1, 2, ..., g. Finally, Shillington suggests

that in the presence of more than one independent variable

problems in grouping may arise, and in this case it may be

wise to perform a different lack of fit test for each

parameter. Following this approach, an example is given

which suggests testing lack of fit for each of the p

independent variables separately may be more powerful than







trying to form groups based on all independent variables at

once.

In summary, all the approaches we have discussed for

testing for lack of fit when replicate observations of the

response are not available at any of the settings of the

independent variables make use of grouping the observed

response values according to similar values of the

independent variables. The observations falling in such

groups are referred to as "pseudoreplicates" or "near

neighbor observations." These pseudoreplicates are used to

estimate the true variance of the observations, a2, but a

completely unbiased estimate of 02 cannot be attained unless

true replicate observations are available. In each case,

the power of the lack of fit testing procedure is reduced

because an unbiased estimate of a2 is not attainable. We

now turn to the use of check points for lack of fit testing.

2.4 Testing for Lack of Fit with Check Points

An alternative to the two approaches to lack of fit

testing already discussed is the method which makes use of

check points. We assume a model of the form E(Y) = Xl as

given in (2.1), is fitted in a response surface system, but

that the true model is of the form E(Y) = X1I + X22 as

given in (2.2). The parameters, a8, in the fitted model are

estimated by ordinary least squares techniques, making use

of the values of the response observed at the design

points. After the model is fitted, values of the response

are observed at additional points in the experimental region







called "check points." The observed response values at the

check points are compared to the values which the fitted

model predicts at these same check points. It is important

to note that the observed values of the response at the

check points are not used in fitting the model initially.

Snee (1977) gives four methods of validating regression

models, one of which is the collection of new data to check

predictions from a previously fitted model. In a designed

experiment these new data take the form of check points.

Snee suggests that the inclusion of a small number of check

points in any designed experiment is a "worthwhile"

procedure.

Scheffe (1958) proposed a test for lack of fit when the

{3,2} simplex lattice design is used for fitting a second

order canonical polynomial model in three mixture

components. It is desired to use the observed value of the

response at (1/3, 1/3, 1/3) as a check point blend. The

test statistic proposed is the t statistic of the form


Y Y
t = (2.6)
[var(Y -)]


where Y is the observed value of the response at the check

point, and Y is the value of the response predicted at the

same point by the second order model which is fitted by

ordinary least squares techniques to the observed response

values at the six design points of the {3,2} simplex

lattice. The response value observed at the point







(1/3, 1/3, 1/3) is not used in fitting the model. Lack of

fit is inferred if the absolute value of the calculated t

value in equation (2.6) is larger than the corresponding

tabled t value.

In the denominator of the t test of equation (2.6), the

variance of the difference Y Y is shown to be



var(Y Y) = var(Y) + var(Y)

= (44/27r)o



when r replicates are taken at each design point. The

estimate of the variance of Y Y is (44/27r)o2, where o2 is

calculated from the replicated response values at the design

points.

Scheffe (1958) also alludes to a test for lack of fit

when several check points are used simultaneously. When

there are k check points, the test for lack of fit is an F

statistic of the form


-1
F = 2 (2.7)
ka

where d' = (Y Y 2 2' ... Yk Yk) and V = o2V

var(d). Formulas are given for the elements of VO in the

special case when the check points are the design points of

the {3,2} simplex lattice. Lack of fit is suspected if the

calculated value of the F statistic given in (2.7) is larger

than the corresponding tabled F value.







Gorman and Hinman (1962) suggest the same t test in

equation (2.6) that Scheffe (1958) suggested for a check

point taken at (1/3, 1/3, 1/3) to test for lack of fit in a

second order polynomial model fitted from a {3,2} simplex

lattice design. They suggest using (1/3, 1/3, 1/3) as the

location of the check point because the observation at this

point may later be used to fit the next more complex model,

the special cubic, if the second order model is found to be

inadequate. They state that in general for the second order

polynomial model as well as higher order models, check

points should be taken in regions of particular interest, of

which there are usually many in any blending study.

Further, they suggest that the number of check points

depends on individual experimental situations--technical

background, precision required, cost of materials and

analyses, and probability of requiring a more complex

model. However, no specific criterion is given by Gorman

and Hinman for selecting the location of the check points.

Gorman and Hinman (1962) indicate that a t test at a

check point other than at (1/3, 1/3, 1/3) takes the same

form as the statistic of equation (2.6),


Y Y
t 1/2
[var(Y) + var(Y)]/2



with the additional condition that if several check points

are taken, say for example k points, the method of checking







the fit is to compute the t value at each location and refer

these calculated t values to the 100(a/2k) percentage point

of the central t distribution rather than the 100(a/2)

percentage point.

Kurotori (1966) gives an example of a mixture

experiment where the response is the modulus of elasticity

of a rocket fuel, which is a mixture of three components,

binder (xl), oxidizer (x2), and fuel (x3). The factor space

of feasible mixtures is a subspace inside the two-

dimensional simplex or triangle where all three components

are present simultaneously. "Pseudocomponents" are defined

and in the pseudocomponent system a special cubic model is

fitted to data collected at the points of the q = 3 simplex

centroid design (Figure 4). A check for adequacy of fit is

made by using three check points and the response values at

the check points are used only for testing the fit of the

model and not for fitting the model initially.

The reason for the choice of the particular check point

locations by Kurotori is that, as he states, "They are the

most remote mixtures from the seven design points." The

lack of fit test is an F statistic of the form


2
F = (2.8)


2 3
where s (Y Y ,) for the i = 1, 2, 3 check points
2 i=l
and 02 is an estimate of measurement error from a previous

analysis. Kurotori admits that the use of the F statistic











x :1
I
(1,0,0)
9 Design Points

O 0 Check Points




SE 3 ') o
S 2 I 112



e(o,,o) a the ser)
2 '2 '2)




Figure 4. Kurotori's rocket fuel example,
xI', x2', and x3' represent pseudocomponents.





in Eq. (2.8) for lack of fit testing may be risky because

the predicted values at the check points are correlated

(correlation of .5), although the observed values are not

correlated. Kurotori suggests individual t tests as

proposed by Scheffe (1958) might be the preferred procedure.

Snee (1971) repeats Kurotori's rocket fuel example

using the same F test for lack of fit as Kurotori and makes

the comment that the Yi's at the check points are

correlated. In stating that the F test is not an exact

test, he nevertheless offers no solution in the form of an

exact test.





39

In summary, only Scheffe refers to an exact F test when

several check points are considered simultaneously for

testing for possible lack of fit of a model fitted in a

mixture space, and his development is limited to the special

case where the check points are the design points used to

fit the model initially. No criterion is proposed by

Scheffe for selecting other locations for the check points.














CHAPTER THREE
AN OPTIMAL CHECK POINT METHOD FOR TESTING
LACK OF FIT IN A MIXTURE MODEL

3.1 Introduction

In Chapter Three we investigate the problem of testing

for lack of fit of a linear model fitted in a mixture

space. The testing is to be accomplished with the use of

check points. We assume that an experimental design is

specified, and that the fitted model is of the form



E(Y) = X-1 (3.1)



where Y is an Nxl vector of observable response values, X is

an Nxp matrix of known constants and rank p, and 81 is a

vector of p unknown regression coefficients. The true model

is assumed to be of the form



E(Y) = X1 + X 2 (3.2)
-1 2-2


where X2 is an Nxp2 matrix of known constants and g2 is a

vector of p2 unknown regression coefficients. Throughout

our development, we will assume that the random vector Y has

the normal distribution with variance-covariance matrix

equal to a IN.








In our investigation we wish to determine the proper

testing procedure to follow in deciding whether the fitted

model exhibits lack of fit. In order to optimize the lack

of fit testing procedure, we will determine the location of

the check points so that the power of the test is maximized.

3.2 Testing for Lack of Fit in the Presence of
an External Estimate of Experimental Error Variation

3.2.1 The Test Statistic

We wish to test the performance or fit of a fitted

model in a mixture space when the true model possibly

contains terms in addition to those in the fitted model.

The fit of the model is to be tested by a test which makes

use of the response values observed at certain locations

called "check points" in the experimental region, by

comparing them to the values which the fitted model predicts

at the same check points. The observed values at the check

points are not used for estimating the coefficients in the

fitted model and are assumed to represent the values of the

true surface at the check points.

Let us define the vector of differences



d = (Y* Y*)



(Y* Y- Y Y* Y* y*)i
1 1' 2 21 k k


where Y* = 1, 2, ..., k are observed response values at

k check points and YV, i = 1, 2, ..., k are response values

predicted at the k check points by the fitted model,







VY = xt'b where b is the ordinary least squares estimator

of al and where x*' is the ith row of X*, the kxp matrix
1 -1
whose columns are of the same form as the columns of X but

with its rows evaluated at the k check points. Note that

if 3 = 0, then E(d) = 0 and if a 2 0, then

2 2 2
E(d) = (X2 X*(X'X)-X'X2)-. Let V represent the

variance-covariance matrix of the random vector d.

Then V = o V where


-1
V = Ik + X*(X'X) X*'



and where Ik is the identity matrix of order kxk.

We assume that an unbiased estimate of a2 is available
^2
and we denote this estimate by a where the subscript ext
Sext'
^2
stands for external, and a is independent of the model

being fitted. The test statistic for the hypothesis of zero

lack of fit H : E(d) = 0 is

-1
d'V d/k
0-
F = (3.3)
^2
ext


(see Scheffe, 1958, p.358). It will be shown later in this

section that the F ratio in Eq. (3.3) possesses either a

central F distribution or a noncentral F distribution,

depending upon whether the true model is represented by Eq.

(3.1) or Eq. (3.2).
^2
The variance estimate a that appears in equation
ext
(3.3) is ordinarily generated from replicated observations








at some of the design points in the experiment. We assume
^2
that aext is a constant multiple of a central chi-square

random variable with v degrees of freedom. This is written

as


^2
a e = SSE /v
ext pure


2 2
= (a 2/v)(SSEpe/
pure


2 2
where SSEpure/ ~ X. Note that SSEpure denotes the

portion of the residual sum of squares due to replication

variation from the fitted model. The residual sum of

squares from the fitted model may be partitioned into

SSEpure and SSLOF only if replicated observations are

collected at one or more design points. For the case where

replicate observations are collected at all of the design

points


n n.
SSE = Z (Yij Yi. )
i=l j=1


where n is the number of distinct design points, n. > 2 is
1
the number of replicates at the ith design point, Yij is the

jth observation at the ith design point, and Y. is the
1.
average of the ni observations at the ith design point.
n
Here SSEpure has v = Z (ni 1) degrees of freedom.
i=l1
When the fitted model and the true model are of the

same form as defined by Eq. (3.1), the quantity d'Vld/a2
(3.lr th quntit d'







possesses a central chi-square distribution (Searle, 1971,

p.57, Theorem 2). However, when the true model is of the
-1 2
form specified by Eq. (3.2), d'V d/a possesses a

noncentral chi-square distribution. Thus when the true

model is of the form in Eq. (3.1),


-1 2 2
d'V- d/ao 2 X
0 Xk'


but when the true model is of the form in Eq. (3.2),


-1 2 2
d'V d/o2 ~ x ,



where in the second case the noncentrality parameter X has

the form



S= E(d)'V E(d)/2a2
1 0


-1 2
= X* X*A)'Vo (X X*A) 2/2
2 2 0 2 2


-1
The matrix A = (X'X) X'X is called the alias matrix and is
2
of order pxp2. In X1, the matrix X* is of order kxp2 and

has the same relationship to X2 as X* has to X.
2
Since SSE /a is statistically independent of
pure
-1 2
d'V d/a then under model (3.1) the test statistic
0 0








-1 2
d'V d/ko
F = 2
SSE pure/v
pure


-1
d'V d/k
= ^2
ext


will have a central F distribution. When the true model

contains terms in addition to those in the fitted model then

F will have a noncentral F distribution. We write these two

cases as



F ~ Fk,
k,v



under model (3.1), and



F ~ F'
k,v;Xi



under model (3.2), where the noncentrality parameter is



x = (X X*A)'V0 (X- X*A)_2/2 2.



3.2.2 The Testing Procedure and an Expression for the Power
of the Test

Given that the form of the fitted model is defined as

Eq. (3.1), the expected value of the numerator of the F

statistic in Eq. (3.3) will depend on the form of the true

model. For the case where the true model is expressed as







Eq. (3.2),



-1
E(numerator) = E(d'V d/k)
0-


2 2
= ( /k)EXx2



= (a /k)(k + 211)


2 2
= 2 + 20 X1/k



= o2 + BA1 /k, (3.4)


-1
where A1 = (X* X*A)'V0 (X X*A). However, when the true

model is Eq. (3.1), 8 = 0 and in this case A = 0 so that
-2 1
2 ^2
E(numerator) = 02. Also a is an unbiased estimator of
ext
o2 and



^2 2
E(a ) = a (3.5)
ext


Therefore the ratio E(numerator)/E(denominator) where
^2
the denominator is ext will equal unity under model (3.1),

that is, when there is no lack of fit. Under model (3.2),

the ratio will be greater than or equal to unity so lack of

fit should be suspected if the calculated F ratio in

equation (3.3) is large. We can thus use an upper tailed

rejection region to reject the hypothesis of zero lack of

fit. The power of the test is










PI k ,v ;X > a; k, v 1
P{Fv > F,; };k,I



where F is the upper 100a percentage point of the
a ;k,v
central F distribution with k numerator degrees of freedom

and v denominator degrees of freedom.

It is worth noting that from Eq. (3.4) and Eq. (3.5)

testing the hypothesis that -2 = 0 is equivalent to testing

the hypothesis that X1 = 0, assuming A1 is positive

definite. Thus testing a null hypothesis of zero lack of

fit using the proposed testing procedure involving the F

ratio in (3.3) may be expressed as a test of the hypotheses:



H: 1 = 0



H: a1 > 0.
a 1


3.2.3 A Method for Locating Optimal Check Points

Once a design for fitting model (3.1) in a mixture

space is chosen and the number of simultaneous check points

is decided on, say k > 1, the next step is to determine

where in the mixture space we should place the k check

points so as to maximize the power of the test for lack of

fit. The location of the check points is to be made

independently of the value of .
-2







The power of the upper tailed F test for lack of fit is

an increasing function of X1 (see Appendix 1 for proof, with

X2 = 0). Therefore, to maximize the power of the test we

maximize the value of X1 defined as







-1
Xl = S2AI82/2a2



where A = (X X*A)'V0 (X* X*A), by properly selecting

the k check points whose coordinates are defined in X*. To

maximize the value of X1, we shall concentrate on the matrix

A1.

The matrix A1 is a square matrix of order p2xp2 and is

a scalar quantity when p2 = 1. By maximizing the scalar

quantity A1 with respect to the k check points, the power is

maximized no matter what the value of .2 Maximizing the
-2
scalar A1 can be accomplished by using The Controlled Random

Search Procedure given by Price (1977). This procedure is

described in Appendix 2. As a computational aid, A1 can be

expressed as


V + (X* X*A)(X* X*A)'
A = V-2 1 (3.6)
1 V


when p2 = 1, where the symbol IBI denotes the determinant of

the square matrix B. Thus the computations reduce to

evaluating two determinants rather than inverting VO (see

Scheffe, 1959, Appendix V, p.417).








When p2 > 1 and A1 is no longer a scalar, maximizing X1

(and thus maximizing the power of the test) cannot be

accomplished without specifying 02- In this case we make

use of a lower bound for 1l (Graybill, 1969, p.330, Theorem

12.2.14(9)) defined as



lmin-/22 < X1



(where min is the smallest eigenvalue of A1) to be used in

place of \l. Hence an approximate solution to the

maximization of X1 will be achieved by finding the k

simultaneous check points (using Price's procedure) that

maximize min' the smallest eigenvalue of A1. In other

words when p2 > 1, and in order to avoid specifying 2' we

seek to maximize a lower bound value for X1. This

maximization does not depend on the value of 02.

There are cases where the matrix A1 is of less than

full rank (less than rank p2) or equivalently where the

matrix A1 is positive semi-definite so that umin will be

equal to zero no matter which check points are selected.

One such case occurs when k < p2 (when the number of check

points is less than the number of parameters in the true

model which are not in the fitted model) since when k < p2



rank(Al) = rank[V I X* X*A)]


rank(X* -X*A),
2







and so rank(Al) < min(k, p2) because the matrix (X* X*A)

is of order kxP2. Therefore when k < p2, the rank of A1 is

at most k so that A1 is of less than full rank. Since umin

must be equal to zero when Al is positive semi-definite, an

alternative method to that of maximizing umin to select

optimal check points must be found when A1 is positive semi-

definite in order to produce a positive lower bound for X1.

In this pursuit, let us write X1 as



1 = -2AA2/202


= PAP'_2/202



= 8[P1:P2] diag[Al, A2=0][P1:P2] '2/2o2


= P-lA 1P 12/202


where A is a diagonal matrix with elements equal to the

eigenvalues of Al, P is an orthogonal matrix whose columns

are orthonormal eigenvectors of AI, A1 and P1 correspond to

the positive eigenvalues of AI, while A2 = 0 and P2

correspond to zero eigenvalues of AI. Then by Theorem

12.2.14(9) in Graybill (1969) we can write



Sz'z/2o2 < 1 (3.7)
min p


where in is the smallest positive eigenvalue of A,, and
min







z = P'S Thus by Eq. (3.7), an approach to maximizing a
- -1-2
positive lower bound for X1 when A1 is positive semi-

definite is to select check points that maximize the

smallest positive eigenvalue of A1. It must be noted,

however, that this method can only be used when

a2 e n C(Pl), where C(P1) denotes the column space of P1
and n C(PI) denotes the intersection of all such spaces

which can be obtained at all possible check points

locations. This is because, in general, z'z in (3.7)

depends on the location of the check points through its

dependency on Pi. If, however, 2 e nC(P ), then

zz = P1 = P = 22 since 'P2 = 0.

It follows that when 2 e n C(P ), mn z'z/2o
2 1 min- -
+ 2 +
= min 2+-/2o and only mn depends on the location of the

check points.

3.3 Testing for Lack of Fit When MSE Is Used
to Estimate Experimental Error Variation

3.3.1 The Test Statistic

In this section we shall show that when an external

estimate of a2 is not available and the residual mean square

(MSE) from the fitted model of the form (3.1) must be used

as an estimate of a2, the test statistic




-l
d'V d/k
F ME (3.8)

MSEpossesses a central F distribution when the true model is

possesses a central F distribution when the true model is







Eq. (3.1), but possesses a doubly noncentral F distribution

when the true model is Eq. (3.2).

In the initial section of this chapter, the quantity

d'V- d/2 was said to possess a central chi-square

distribution or to possess a noncentral chi-square

distribution, depending on whether the true model was

specified by Eq. (3.1) or Eq. (3.2). Now, the residual sum

of squares from the fitted model is defined as


N
SSE = E (Y Y)
i=l

-1
= Y'(I X(X'X) X')Y



and it is easy to show (Searle, 1971, p.57, Theorem 2) that

SSE/a2 possesses a central chi-square distribution if the

true model is Eq. (3.1), but under model (3.2), SSE/a2

possesses a noncentral chi-square distribution. This is

expressed as


2 2
SSE/a xp
N-p


under model (3.1), and



2 2
SSE/a X.
N-p,12







under model (3.2), where the noncentrality parameter \2 is


2
X2 = -(X2 XA)'(X2 XA)2/2o0



The distributional form of the test statistic in Eq.

(3.8) is derived by knowing that the quantities

d'Vo d/o2 and SSE/o are statistically independent (see

Appendix 3), so that


d'vold/ko2



-12
MSE/o



d'V d/k
MSE


is distributed as a central F when the true model is Eq.

(3.1), but when the true model is Eq. (3.2) the F ratio is a

doubly noncentral F, that is, under model (3.1),



F ~ F
k,N-p


and under model (3.2),



F~ F"
k,N-p;X ,X2


3.3.2 The Rejection Region and its Relation to the Power of
the Test

In Appendix 1 it is shown that if k, N-p, and X2 are

fixed, then the power of the F test using the ratio (3.8) is







a function of the location of the rejection region (upper

tailed or lower tailed) of the test. The power increases

with increasing values of the numerator noncentrality

parameter, X1, when the test is an upper tailed test. The

power decreases with increasing values of X1 when the test

is a lower tailed test. This means that to study ways of

increasing the power of the test, we have to determine

whether the test is an upper tailed test or a lower tailed

test. Similarly, for fixed values of k, N-p, and X1, the

power of the F test is a decreasing function of X2 for an

upper tailed test, and is an increasing function of X2 when

the F test is a lower tailed test (Scheffe, 1959, p. 136-

137).

To decide if the test is an upper tailed test or a

lower tailed test, we recall from Section 3.2.2 that if the

true model is Eq. (3.1) then the expected value of the

numerator of the F statistic in (3.8) can be written as



E(numerator) = 02,



and if the true model is Eq. (3.2),




2 2
E(numerator) = a + 2a X1/k (3.9)



= 02 + QA 2/k








-1
where the P2xp2 matrix A is A = (X* X*A)'V (X* X*A).

Similarly, it can be shown that if the true model is Eq.

(3.1), the expected value of the denominator of the F

statistic in (3.8), where the denominator equals MSE, is



E(denominator) = E(MSE)



= 2,



but if the true model is Eq.(3.2),



E(denominator) = E(MSE)



[o2/(N p)]Eyp,2


= [a2/(N p)][N p + 2X2]



= a2 + 22 2X/(N p) (3.10)



= 2 + a2A2 2/(N P)



where the P2xP2 matrix A2 is A2 = (X2 XA)'(X2 XA). Thus

the ratio E(numerator)/E(denominator) will equal unity if

the true model is Eq. (3.1), but if the true model is Eq.

(3.2), the ratio is greater than unity if B_~A12/k >

8_A 2 2/(N p). In this latter case we reject the null







hypothesis of zero lack of fit if the calculated value of

the F ratio in (3.8) is large. An upper tailed rejection

region seems reasonable for this test. When the true model

is Eq. (3.2), and if a2Al2/k < 2A2 2/(N p), then a lower

tailed rejection region is preferred.

3.3.3. A Method for Locating Optimal Check Points

Given a design for fitting a model of the form in Eq.

(3.1) in a mixture space (note that fixing the design fixes

A2 and (N p)), and given the number of simultaneous check

points desired, k > 1, we now wish to determine where in the

mixture space the k check points should be located so as to

maximize the power of the F test for lack of fit, where the

test statistic is given in Eq. (3.8). We also wish to

position the optimal check points in a manner that is

independent of the values of the elements in 2 .
-2
The case of an upper tailed test. To help us find k

simultaneous check points that maximize the power of an

upper tailed test, we shall make use of the fact that the

power is an increasing function of Xi. Therefore to

maximize the power of the upper tailed F test, we shall seek

the locations of the k check points that maximize X1.

As in the case considered in Section 3.2.3, where the

test statistic had a noncentral F distribution, if the

number of extra terms in the true model is p2 = 1, then

maximizing X1 is equivalent to maximizing the scalar A1.

However, as before, if p2 > 1, then the P2xP2 matrix A1 is

not a scalar and we will have to approximate the







maximization of X1 by maximizing a lower bound for X1. This

is done by finding the maximum value of min' the smallest

eigenvalue of A1, since



minm-2-/22 X1


When the number of check points is less than the order

of the square matrix A1, that is, k < p2, then rank(A1) <

min(k, p2), and A1 will have Umin = 0. For this case, we

again try to maximize the smallest positive eigenvalue of A1

which we denote by in' while remembering from Section
min'
3.2.3 that this technique is limited to situations where

B2 e nC(P1)
-21
The case of a lower tailed test. To find k check

points to maximize the power of a lower tailed test, we make

use of the fact that the power of the lower tailed F test

increases as X1 decreases. Then if P2 = 1 and A1 is a

scalar quantity, X1 can be minimized with respect to the k

check points by finding the check points that minimize A1.

If P2 > 1, then by Theorem 12.2.14(9) in Graybill (1969), we

see that an upper bound for X1 is



X1 max_22 /202 (3.11)



where umax is the largest eigenvalue of Al. An approximate

solution to minimizing X1 in (3.11) can be achieved by

minimizing max. It is not necessary to treat the case
-'max







of k < p2 separately here, although X1 will equal zero if

_2 is in the column space of P2, where P2 is the matrix
whose columns are orthonormalized eigenvectors corresponding

to the zero eigenvalues of the matrix A .

3.3.4 Determining Whether the Test Is Upper Tailed or Lower
Tailed

The procedures outlined in Section 3.3.3 produce a set

of k check points that simultaneously maximize the power of

the upper tailed test as well as a second set of k check

points that simultaneously maximize the power of the lower

tailed test. The check points that are selected maximize

the power, given A2, k, and N p without specification

of a2' except that when A1 is positive semi-definite we
require that 82 e n C(P ).

It is now necessary to decide which of our two

candidates will be used for a lack of fit test. To choose

between the upper tailed test and the lower tailed test, let

us consider the quantity



R = [A /k] [A2/(N p)].



If R is positive definite when the true model is Eq. (3.2),

then no matter what the value of 82 is, the ratio

E(numerator)/E(denominator) will be greater than unity,

implying an upper tailed test is to be used. Similarly, if

R is negative definite, then a lower tailed test should be

used. Finally, if R is not definite, then neither an upper

nor a lower tailed test is implicated and further







investigation is necessary. The criterion of R = [A /k] -

[A2/(N p)] may yield any of the four following cases.

Case 1. If R = [A /k] [A2/(N p)] is positive

definite when A1 is generated by the k optimal upper tailed

test check points, and R is not negative definite when A1 is

generated by the k optimal lower tailed test check points,

then we recommend that the check points be used that yield

the optimal upper tailed test with an upper tailed rejection

region.

For Case 1 it is necessary for A1 to be positive

definite (see Appendix 4). Since A1 is a square matrix of

order p2xP2 with rank(A ) < min(k, p2), then A1 can be

positive definite only if k > p2. Thus, there must be at

least p2 check points for Case 1 to hold, where p2 is the

number of terms in the model of Eq. (3.2) that are not in

the model of Eq. (3.1).

From inspection of equations (3.9) and (3.10), it is

apparent that the testing for lack of fit in Case 1 is

equivalent to testing the hypothesis




1 2
H0 0 N p = (3.12)


against the alternative




1 2
H > 0
a k N p





60

since R = [A /k] [A2/(N-p)] is positive definite when the

true model is Eq. (3.2). In Appendix 5(a) it is shown that

under Case 1, the hypothesis given by (3.12) is equivalent

to the hypothesis



H X = X = 0.



Case 2. In Case 2 we assume that R = [A /k] -

[A2/(N p)] is not positive definite for the k optimal

upper tailed test check points, but that R is negative

definite for the k optimal lower tailed test check points.

Here we recommend that the lower tailed test check points be

used with a lower tailed rejection region.

It is necessary for A2 to be positive definite for Case

2 to occur (see Appendix 4). However, A1 need not be

positive definite, and so k need not be greater than p2. In

Case 2 then, it is possible that lack of fit may be tested

with only one check point.

By inspection of equations (3.9) and (3.10), a

hypothesis of no lack of fit is equivalent to


X1 X
1 2
Hk = 0 (3.13)
0' k N p


while the alternative hypothesis that lack of fit is present

is equivalent to


X1 X
1 2
H N < 0
a k N p







since R = [A /k] [A2/(N p)] is negative definite. In

Appendix 5(b) it is shown that the hypothesis given by

(3.13) is equivalent to the hypothesis



H0: = 2 = 0.


Case 3. We assume R is positive definite for the k

optimal upper tailed test check points, and R is negative

definite for the k optimal lower tailed test check points.

Hence either an upper or lower tailed test may be considered

as a possible test for lack of fit. If the quantity
2
_'_2//o can be specified, then the minimum power for both
the optimal upper and optimal lower tailed tests can be

approximated, and the test with the greater minimum power is

recommended. In Appendix 4 it is shown that Case 3 can

occur only when A1 is positive definite for the upper tailed

test. Thus Case 3 can only occur when there are at least p2

check points.

The minimum power of the upper tailed test may be found

by calculating



P IF" > F ), (3.14)
k,N-p;1 IL' 2U a;k,N-p


where F;k,N-p is the upper 100a percentage point of the

central F distribution,






2
IL = /min-2/2a2


and



X2U = max2-2/2 02



where min is the smallest eigenvalue of A1 and max is the
mmn 1 max
largest eigenvalue of A2. Formula (3.14) yields a

conservative lower bound for the power of the optimal upper

tailed test. Note that A1 is generated using the optimal

upper tailed test check points. The cumulative distribution

function of F" can be approximated by multiplying the

cumulative probabilities of the central F distribution by a

constant (Johnson and Kotz, 1970, p.197). This

approximation is described in Appendix 6. Other

approximations for F" (such as the Edgeworth series

approximation suggested by Mudholkar, Chaubey, and Lin,

1976) exist which are generally more accurate, but we chose

to use the approximation given in Johnson and Kotz (1970,

p.197) due to its simplicity. Additionally, the

approximation of Mudholkar, Chaubey, and Lin (1976) produced

negative probabilties when only one degree of freedom was

available in either the numerator or denominator of F".

This problem was avoided by using the approximation given by

Johnson and Kotz (1970).

The minimum power of the optimal lower tailed test can

be approximated similarly (if B~22/o2 is specified) by







calculating


P IFF"


< (l-a);k,N-pp


where


X lu = maxJ.2 2/2o 2


and


S2L = 6min-22/20 2


with Umax equal to the largest eigenvalue of A1 and 6min
max n
equal to the smallest eigenvalue of A2. Note that A1 is

generated by using the optimal lower tailed test check

points. For the lower tailed test, A1 may be positive semi-

definite, and if 82 is in the column space of P2 then A1 = 0.

In Case 3, the upper tailed test is a test of



HO: 1 = = 0




0 1 2
H 2 > 0
a k N- p


while the lower tailed test is a test of







H: X =X2 0
H0 1 2 = 2



1 2
H -P < 0.
a k N P < 0


Case 4. In Case 4 we assume that R = [Al/k] -

[A2/(N p)] is not positive definite for the k optimal

upper tailed test check points and R is not negative

definite for the k optimal lower tailed test check points.

Here it is useful to write the difference between the

expected value of the numerator and the expected value of

the denominator of the F ratio in (3.8) as



al[A /k A2/(N P)] = -sns'S
-2 1 2 -2 = -2


= 8[S1:S2:S3] diag[l2'r2=0,03 [Sl:S2:S3] '2



= 8_2S"l I2 + 2S3 332


where 0 = diag(ll, 02' 23) is a diagonal matrix consisting

of the eigenvalues of R, 01 is a diagonal matrix of the

positive eigenvalues of R, 02 is a diagonal matrix of the

zero eigenvalues of R, and 03 is a diagonal matrix of the

negative eigenvalues of R. The orthogonal matrix S can be

expressed as S = [S1:S2:S3], where the matrices SI, S2, and

S3 have columns which are orthonormalized eigenvectors

corresponding to nl, 02, and 03, respectively.







In Case 4, neither the optimal upper tailed test nor

the optimal lower tailed test is applicable for all values

of _2 For completeness, we note that Case 4 actually

consists of nine subcases, where R may be positive semi-

definite, negative semi-definite, or indefinite for either

the optimal upper tailed test or lower tailed test check

points. These subcases are listed in Table 2.




Table 2. Nine Subcases of Case 4.


R--Upper R--Lower
Subcase Tailed Test Tailed Test

1 PSD PSD
2 PSD NSD
3 PSD I
4 NSD PSD
5 NSD NSD
6 NSD I
7 I PSD
8 I NSD
9 I I

PSD = positive semi-definite, NSD = negative semi-
definite, I = indefinite.





If _2 lies in the column space of S2, then 8'[A /k -

A2/(N p)]s82 is zero, and therefore lack of fit is not

testable with either an upper or lower tailed test. A

sufficient condition for the test for lack of fit to be

upper tailed in Case 4 is that 0 be in the column space
-2
of [S1:S2], but not entirely in the column space of S2. In

this case








;[A1/k A2/(N p)]_2 = 2S01812 + 2S303S32


= 2S1~lS 2 + 0


= S~Sl2lSI2,


and 8_[ Al/k -
indicating an
condition for
that _2 be in
in the column


A2/(N )]_2 will be greater than zero,
upper tailed test. Similarly, a sufficient
the test for lack of fit to be lower tailed is
the column space of [S2:S3], but not entirely
space of S2. Then


2[A1/k A2/(N P)]g2 = 0 + 2S33S3_2


= 2S3 3S 2


which makes _2[Al/k A2/(N p)] 2 less than zero,
indicating a lower tailed test.

To determine whether 2 is in the column space of

[S :S2], let us define the augmented matrix

Q1 = [f2:Sl:S2] If Q!Q1 has a zero eigenvalue, then 02 is
in the column space of [S :S2]. Similarly, if we define

Q2 = [82:S2] and Q3 = [82:S2:S3]' then a2 is in the column
space of S2 if Q'Q2 has a zero eigenvalue, and 2 is in the
column space of [S2:S3] if Q3Q3 has a zero eigenvalue.







Given that we are in a particular subcase of the nine

subcases described in Table 2, we recommend that lack of fit

be tested with the upper tailed test check points if it is

determined that 2 is such that '[A /k A2/(N -p)]2 is

positive when A1 is generated from the upper tailed test

check points. Likewise, for the same given subcase, if the

value of 02 of interest is determined to produce a negative

value for i2[Al/k A2/(N P)]@2 when A1 is generated from

the lower tailed test check points, then we recommend that

lack of fit be tested with the lower tailed test.

We see then that Case 4 is an undesirable situation in

practice, since, in order to test for lack of fit, we must

assume a priori that any lack of fit is due to a nonzero

value of 82 that produces an upper tailed or lower tailed

rejection region. However, it would seem rare that such

knowledge would be available.

3.4 Examples

We now present several examples to illustrate the

technique for locating optimal check points to be used in

testing for lack of fit in a mixture model.

3.4.1 Theoretical Examples

Example 1. In this example a second order canonical

polynomial model is fitted in three mixture components using

the {3,2} simplex lattice design, which is presented in

Figure 1 of Chapter 1. The true model is assumed to be the

special cubic model containing the term 123x x 23 in

addition to the six terms of the fitted model. The expected







values of the response at the six design points are assumed

to be represented by the fitted model in the form



E(Y) = X01,



but with the true model the expectations are written as



E(Y) = X 1 + X22 ,



where X is a 6x6 matrix with rows that define the

coordinates of the six design points and columns that

correspond to the six terms in the fitted model (xi, x2, x3,

xlx2, xlx3' x2x3)' 8a is the 6x1 vector of regression

coefficients (81, 82, 83, 812, 813, 823), X2 is a 6x1

column vector containing the values of the term xlx2x3 at

the design points, and 82 is the single regression

coefficient 8123'

The {3,2} simplex lattice design consists of only six

design points, and since six parameters are estimated in the

second order fitted model, there are no degrees of freedom

remaining for obtaining an estimate of the experimental

error, 02. We assume therefore that an external estimate of

a is available, a2 which will be used in the denominator

of the lack of fit F statistic given in Eq. (3.3).

Since there is one term in the true model in addition

to those in the fitted model, that is p2 = 1, we know that

in order to locate k simultaneous check points that maximize







the power of the test for lack of fit it is necessary to

maximize the scalar quantity



-1
A = (X* X*A)'V (X* X*A)
1 2 0 2


with respect to the coordinates of the k check points. Here

X* is a k-element column vector with ith element equal to

the value of x* x* x* at the ith check point, X* is a kx6
il i2 i3
matrix with ith row equal to the value of (x* x* x*
il 12' i3'
S x* x* x xt ) at the ith check point,
11 1i2' 11 Xi3' 12 i3
A = (X'X)-X'X2 is the 6x1 alias vector, and

V = Ik + X*(X'X)- X*'. This maximization is accomplished

by use of the Controlled Random Search Procedure (Price,

1977), which is described in Appendix 2.

When only a single (k = 1) optimal check point is

desired the Controlled Random Search Procedure locates a

point (x*, x*) which maximizes


-1
A1 = (X* X*A)'V0 (X2 -X*A),



where


S= x*x*x* = *x*( x* -
2 123 12 1 2







X* = (x*, x*, x* x* x*x*, x*x*)



= (x*, X* (1 x* x*), x*x*, x*(l *),
1 2 1 2 12' 1l 1 2


x*(l x* x*)),
2 1 2

-l
and V0 = 1 + X*(X'X)- *'. The value of A1 is calculated

using the formula of Eq. (3.6). Following this procedure,

we find that the single check point that maximizes A1, and

thus maximizes the power of the test, is the centroid of the

triangular factor space (1/3, 1/3, 1/3). The value of A1 at

this centroid point is A1 = 0.00084.

When the Controlled Random Search Procedure is used to

locate k = 2 simultaneous check points that maximize A,, the

centroid (1/3, 1/3, 1/3) is selected twice, and A1 =

0.00121. For three simultaneous optimal check points, the

centroid is selected three times, and A1 = 0.00142.

To test whether the second order model exhibits lack of

fit, when we suspect the special cubic model is the true

model, we form the F ratio



-1
d'V d/k
-0-
F = ^ 2
^2
ext




with the single check point (1/3, 1/3, 1/3) where d =

Y* Y*, Y* is the observed response, Y* is the response
1 1 1 1







predicted by the second order fitted model at (1/3, 1/3,

1/3), and V0 = 1 + (1/3, 1/3, 1/3, 1/9, 1/9, 1/9)(X'X)-1

(1/3, 1/3, 1/3, 1/9, 1/9, 1/9)'. If the calculated value of

the F ratio exceeds F where v equals the number of

degrees of freedom associated with a then we reject the
ext
null hypothesis that the second order model is the true

model in favor of the alternative hypothesis that the

special cubic model is the true model. Equivalently, we

reject H: I1 = 0 in favor of Ha: 1 > 0. For k = 2 or

k = 3 check points, the value of the F ratio is calculated

using the observed and predicted responses at the two or

three replicates at the centroid. The hypothesis

H0: X = 0 is rejected in favor of Ha: X > 0 if F

exceeds F;kv
a ;k,v
Example 2. In Example 2 we illustrate the second of

the four cases that could arise when MSE is used as an

estimate of 02 in the lack of fit test statistic (see

Section 3.3.4). We again fit a second order canonical

polynomial model in three mixture components, and assume the

true model is special cubic. The design to be used is the

q = 3 simplex centroid design, which consists of seven

design points, and is illustrated in Figure 2 of Chapter 1.

There are six parameters to be estimated and seven

design points hence one degree of freedom can be used to

calculate MSE. We shall use MSE to estimate a2. Optimal

upper and lower tailed test check points must be located,

and then a decision is made as to which test should







be used. The actual testing for lack of fit involves the F

statistic in (3.8).

As in Example 1, p2 = 1, since there is one term in the

true model in addition to those in the fitted model. Thus

Al is a scalar whose value we seek to optimize with respect

to the desired number of check points, k. When only a

single check point is sought for the purpose of testing lack

of fit, the Controlled Random Search Procedure has two

functions. First, the procedure is used to locate the

optimal candidate check point for an upper tailed test by

locating the check point that maximizes the scalar Al.

Secondly, the procedure is used to locate the optimal

candidate check point for a lower tailed test, which is

accomplished by locating the point that minimizes A1. The

quantity R = [Al/k] [A2/(N p)] is then calculated to

determine whether the upper or lower tailed test will be

used. If R is positive for the candidate check point for an

upper tailed test, then the test is upper tailed, and the

test is lower tailed if the candidate check point for a

lower tailed test produces a negative value for R. Note

that A2 = (X2 XA)'(X2 XA) is fixed once the design is

specified, since A2 does not depend on the check points.

Using the Controlled Random Search Procedure it is found

that the maximum value of A1 occurs at (xl, X*, x*) = (1/3,

1/3, 1/3), which will be the location for the check point

for the upper tailed test. Calculating A1 at this centroid







point, we find that R = [A /k] [A2/(N p)] = [(3.7258

x 10-4)/l] [(8.4175 x 10-4)/l] = -4.6917 x 10-4. Since R

is negative, the test is not upper tailed.

Using the Controlled Random Search Procedure to

minimize AI, we find that a subregion of the factor space

exists in which all points yield a near minimum value for

A1. We choose the point (0.0189, 0.9269, 0.0542) at random

from this subregion to be used as the optimal candidate for

a lower tailed test. Here R = 0 [(8.4175 x 10-4)/] =

-8.4175 x 10-4.

Since R is negative for both the optimal upper tailed

test check point and for the optimal lower tailed test check

point, we have Case 2 of Section 3.3.4. The upper tailed

test check point is disregarded, and the lower tailed test

check point (0.0189, 0.9269, 0.0542) is used to test for

lack of fit. If the calculated F ratio,

-1
d'V d
F =
MSE


is less than F ( );,then H: X = X = 0 is rejected in

favor of Ha: [X1/1] [X2/1] < 0, that is we conclude that

the second order model exhibits lack of fit, and the true

model is special cubic.

When two simultaneous check points are desired for

testing lack of fit, we can again use the Controlled Random

Search Procedure to locate the optimal settings. To

maximize the scalar A,, we find that both check points







should be selected at (1/3, 1/3, 1/3), for an upper tailed

test. With our calculations R = [(5.8275 x 10-4)/2] -

[(8.4175 x 10-4)/1] = -5.5038 x 10-4, but since R is

negative, the test is not upper tailed.

Minimizing A1 to locate two optimal lower tailed test

check points yields a subregion in the factor space of

optimal check points. The pair of check points (0.3749,

0.5752, 0.0499) and (0.5332, 0.4169, 0.0499) is selected at

random from this subregion, and these check points yield

R = 0 [(8.4175 x 10-4)/1] = -8.4175 x 10-4.

Since R is negative for the upper tailed test points

and the lower tailed test points, we have Case 2 of Section

3.3.4 again and the lower tailed test check points are used

to test for lack of fit. The hypothesis H0: X1 = 2 = 0 is

rejected in favor of Ha : [1/2] [x2/1] < 0 if the cal-
-1
culated value of F = (d'Vo d/2)/MSE is less than

F(1-a);2,1, in which case we say lack of fit of the model is

present.
^2
If an external estimate aext had been available for
ext
this example, then the optimal upper tailed test check

points could have been used in the F ratio,
-1 ^2
F = (d'V d/k)/o t and lack of fit would then be detected
0 k exta

if the calculated value of F exceeded F ,
a;k,v
Example 3. Example 3 illustrates the procedure for

locating optimal check points when there are two terms in

the true model in addition to those in the fitted model. A

second order canonical polynomial model in three mixture







components is fitted using a q = 3 simplex centroid design.

The true model is assumed to contain eight terms, six of

which are the same terms as in the fitted model, with the

additional two terms being the third order terms

612x x2(X1 x2) and B123x x2x3. As in Example 2, there is

one degree of freedom for MSE which is used to estimate
2 T-
2. The test statistic, F = (d'V d/k)/MSE, is given in

equation (3.8).

Since p2 = 2 and A1 is a 2x2 matrix, locating the

optimal upper tailed test check points by the procedure of

maximizing X1 is assisted by the maximizing of a lower bound

for 1 namely maximizing u min /2a where u is the
I min-2-2 mm
smallest eigenvalue of A1. Since 82 and a2 are unknown,

this is equivalent to maximizing min For min to exceed

zero, it is necessary that A1 be of full rank, and since

rank(A1) ( min(k, p2), it is necessary to select k > 2 check

points. If A1 is less than full rank, and thus is positive

semi-definite, only a subset of possible values of 2 could

be considered to make it possible to test for lack of fit

with an upper tailed test.

Using the Controlled Random Search Procedure, the

points that maximize min are found to be (0.418, 0.277,
mln
0.305) and (0.277, 0.418, 0.305). These points are thus

optimal candidates for upper tailed test check points. At

these check points we have umin = 5.1623 x 10-4, A1 =

diag[5.1623 x 10-4, 5.1916 x 10-4], A2 = diag[0, 8.4175 x

10-4], and R = [A1/2] [A2/I] = diag[2.5811 x 10-4, -5.8217




76

x 10-4]. Since the eigenvalues of R are -5.8217 x 10-4 and

2.5811 x 10-4, R is indefinite. Following the suggested

procedure for Case 4 of Section 3.3.4, we note that an upper

tailed test for lack of fit exists if the value of

'2 = [512' 8123]' is in the column space of [Sl:S2] but not
entirely in the column space of S2, where S1 is the matrix

whose columns are the orthonormalized eigenvectors of R

corresponding to the positive eigenvalues of R, and S2 is

the matrix whose columns are the orthonormalized

eigenvectors of R corresponding to the zero eigenvalues of

R. Since R has no zero eigenvalues in this example, S2 does

not exist, but S1 is the column vector, S1 = [1,0]'. Thus

if 2 is of the form 2 = [612 0]', where 12 0, then

a2 is in the column space of S1 and the test is upper
tailed.

The matrix A2 has rank one and therefore is positive

semi-definite. Hence it is impossible to locate two check

points that minimize umax and also make R = [A1/2] [A2/1]

negative definite (see Appendix 4), that is, it is

impossible to find a lower tailed test that is capable of

testing lack of fit for all values of 2. However, if we

use the Controlled Random Search Procedure to locate two

check points that minimize an upper bound for X1 which is

umax S2 /2 2, then by minimizing max, we find that any of
the check points in a particular subregion of the factor

space yield a near minimum for umax. One pair of points in

this subregion is selected as the points to be used as








optimal lower tailed test check points, namely the pair

consisting of the point (0.053, 0, 0.947) replicated twice.

Replicating this check point, we find max = 7.3900

x 10-11, A = diag[0, 7.3900 x 10-11], A2 = diag[0, 8.4175 x

10-4], and R = [A1/2] [A2/1] = diag[0, -8.4175 x 10-4].

The eigenvalues of R are 0 and -8.4175 x 10-4 implying that

R is negative semi-definite. The values of 8 that are in
-2
the column space of [S2:S31 but not entirely in the column

space of S2 will provide a lower tailed test. Here, [S2:S3]

= diag[l,l] and S2 = [11,0]'. Thus, the test for lack of fit

is lower tailed if 8123 0.

For values of 8 that produce an upper tailed test we
-2
use the check points (0.418, 0.277, 0.305) and (0.277,

0.418, 0.305) with the F ratio



-1
d'Vo d/2
F =
MSE


and conclude there is lack of fit if the calculated value of

F exceeds Fa;21. For values of a2 that produce a lower

tailed test, we use two replicates of the check point

(0.053, 0, 0.947), and'conclude there is lack of fit if F is

less than F(1-a);2,1, where again F is calculated by
-1
F = (d'V d/2)/MSE.
0O
Example 4. Example 4 illustrates Case 3 of Section

3.3.4 in which MSE is used to estimate 02 in the lack of fit

test statistic. A second order canonical polynomial model







in three mixture components is fitted using the {3,3}

simplex lattice design, which appears in Figure 5. The true

model is assumed to be special cubic, thus p2 = 1 and A1 is

a scalar. The {3,3} design consists of ten design points

and since there are six parameters to be estimated in the

fitted model, 02 can be estimated by MSE with N p = 10 6

= 4 degrees of freedom.

We first suppose that a single check point is to be

used to test for lack of fit. Using the Controlled Random

Search Procedure we find the single check point that

maximizes the scalar


-1
A= (X X*A)'V (X* X*A)
1 2 0 2


is located at the centroid of the simplex factor space.

Thus (x*, x*, x*) = (1/3, 1/3, 1/3) is the optimal candidate

for an upper tailed test check point. At this centroid

point, A1 = 4.9076 x 10-4. For the {3,3} design the scalar

quantity A2 = (X2 XA)'(X2 XA) is fixed and is equal to

A2 = 9.4062 x 10-4 and thus, R = [Al/k] [A2/(N p)] =

[(4.9076 x 10-4)/1] [(9.4062 x 10-4)/4] = 2.5560 x 10-4

The point that is the optimal candidate for a lower

tailed test check point is chosen randomly from a subregion

of points in the factor space, in which all points minimize

A1. The point selected has the value (x*, x*, x*) = (0.560,

0.410, 0.030). Here A1 = 9.6590 x 10-7 and R = [(9.6590 x

10-7)/1] [(9.4062 x 10-4)/4] = -2.3419 x 10-4.












(1,0,0)







3 3 a0o 0'x
( 0)/ 3 3

3' 3'3

(0o,,0) --- (0,0,1)
x2 0, ) 3 3' ) x3= 1

Figure 5. The {3,3} simplex lattice design.





Since R is positive for the optimal upper tailed test

check point (1/3, 1/3, 1/3) and R is negative for the

optimal lower tailed test check point (0.560, 0.410, 0.030)

we are in Case 3 of Section 3.3.4. Either the upper or

lower tailed test could be used to test for lack of fit, but

if the quantity B B2/o2 can be specified, then we will

choose to use the test that has greater minimum power, since

greater power means that we are more likely to detect lack

of fit when in fact lack of fit exists. In this example

-2 123*
For illustrative purposes, we arbitrarily choose
2
8'8 /C2 = 2000, so that an approximate conservative lower
bound for the power of the upper tailed test is found by







calculating



P {F" > F k
k,N-p;XIL' 2U a;k,N-p


where F;k,N-p is the upper 100a percentage point of the

central F distribution, k is the number of check points, N

is the total number of response observations, p is the

number of parameters in the fitted model,
2 2
IlL = UminB-2/2a and 2U 6 axI/20 The

quantity umin is the smallest eigenvalue of AI, where A1 is

evaluated at the optimal upper tailed test check point.

Since A is a scalar, mn = A. Likewise, 6 is the
1 mim 1 fmax
largest eigenvalue of A2, and since in this example A2 is a

scalar, 6max = A2. In this example we have k = 1, N p =

10 6 = 4, AlL = UminS2/2a2 = (4.9076 x 10-4)(2000/2)
-1 2
= 4.9076 x 101, and 2U = 6 m /2a2

= (9.4062 x 10 )(2000/2) = 9.4062 x 101. Using the

approximation to the cumulative probabilities of the doubly

noncentral F distribution given by Johnson and Kotz (1970,

p.197) which is described in Appendix 6, and taking a = .05,

we find that a conservative lower bound for the power of the

optimal upper tailed test is approximately equal to .0649.

The minimum power for the optimal lower tailed test is
2
approximated (assuming 8_Y2/a = 2000) by calculating



P F" < FN
k,N-p;lU' 2L (1-a);k,N-p







2
The quantities A1U and 2L are taken as A U = maxJ2 /2o
lJU 2L lU max-z22
and 2L = 6 mn /2o where max is the largest eigenvalue
2L min-2-2 max
of A1 with A1 calculated using the optimal lower tailed test

check point, and where 6 in is the smallest eigenvalue of

A2. Since A1 and A2 are scalars, max = A and 6 min = A .
2 1 2 max 1 mm 2
In this example, k = 1, N p = 4,

AlU = (9.6590 x 10 )(2000/2) = 9.6590 x 104, and
-4 -i
x2L = (9.4062 x 10 )(2000/2) = 9.4062 x 101 Again if the

approximation to the doubly noncentral F distribution given

in Johnson and Kotz is used, an approximate conservative

lower bound for the power of the optimal lower tailed test

is .0555.

Having specified a 2/02 = 2000, the optimal upper

tailed test is chosen over the optimal lower tailed test,

because the approximate minimum power of the upper tailed

test is greater than the approximate minimum power of the

lower tailed test. Using the optimal upper tailed test

check point (1/3, 1/3, 1/3) in the test statistic

-1
d'V d
F =
MSE


we conclude that lack of fit is significant if the

calculated value of F exceeds F a1,4 in which case we
a;1,4'
reject H0: XI = 2 = 0 in favor of Ha: A/l 12/4 > 0.

When two simultaneous check points are used for testing

lack of fit, the Controlled Random Search Procedure locates

the optimal upper tailed test and optimal lower tailed test







check points. It turns out that two replicates at (1/3,

1/3, 1/3) maximize Al, and are used as optimal check points

for an upper tailed test. The value of R = [A,/2] [A2/4]

is [(7.9210 x 10-4)/2] [(9.4062 x 10-4)/4] = 1.6090 x

10-4.

In searching for two optimal lower tailed test check

points, again a subregion of the factor space is found in

which any of the points nearly minimize A From this

subregion are chosen the points (0.6386, 0.3263, 0.0351) and

(0.7257, 0.2421, 0.0322) resulting in a value of R = [A,/2]

- [A2/4] of [(1.5216 x 10-9)/2] [(9.4062 x 10-4)/4] =

-2.3516 x 10-4.

In conclusion, when two simultaneous check points are

used in the test for lack of fit in this example, R is

positive for the optimal upper tailed test and R is negative

for the optimal lower tailed test, and we have Case 3 of
2
Section 3.3.4. Selecting 002 /o = 2000 arbitrarily, we

found the approximate lower bound for the power of the upper

tailed test to be .0504, and the approximate lower bound for

the power of the lower tailed test to be .0612. Since the

power is higher with the lower tailed test it is our choice

for testing lack of fit when two check points are used

simultaneously. Lack of fit is detected and we reject

H0: 1 2 = 0 in favor of Ha: [1 /2] [ 2/4] < 0 if the F
-1
ratio, F = (d'V0 d/2)/MSE, using the optimal lower tailed

test check points (0.6386, 0.3263, 0.0351) and (0.7257,

0.2421, 0.0322) is calculated to be less than F1-a);2,4'
(i-a);2 ,4"







3.4.2 Numerical Examples

Numerical Example 1. In this example we illustrate

numerically some of the findings in the first theoretical

example of Section 3.4.1. Data that were collected in a

rocket fuel experiment (Kurotori, 1966) will be used to

investigate the power of the lack of fit F test. The test

is set up to detect the inadequacy of a fitted second order

canonical polynomial model when the true model is special

cubic. Calculated values of the power of the test which

detects lack of fit through large values of

-1
d'V0 d/k
F =
^2
ext

will be compared for several check point locations, includ-

ing the location (1/3, 1/3, 1/3) at which the power was

found to be maximum in Example 1 of Section 3.4.1.

In Kurotori's experiment the modulus of elasticity (Y)

of a rocket fuel is expressed as a function of the

proportions of three components--binder (xl), oxidizer (x2),

and fuel (x3). Since lower bounds are placed on the

component proportions xl, x2, and x3, in the form of

0.20 < xl, 0.40 < x2, and 0.20 < x3, pseudocomponents (x!)

are defined in terms of the original components in the form

of x1 = (xl 0.20)/(1 .80), x' = (x2 0.40)/(1 .80),

and x' = (x3 0.20)/(1 .80). The true special cubic

model in the pseudocomponents, which is obtained by using

the data at the seven points of the simplex centroid design







in the pseudocomponent system, is



E(Y) = 2350x' + 2450x' + 2650x' + Ox'x'
1 2 3 1 2


+ lOOOx'x3 + 1600x'x' + 6150x'x'x'.



The second order canonical polynomial model that is fitted

to the six boundary points only, and which will be tested

for lack of fit, is given by



Y = 2350x' + 2450x' + 2650x'
3


+ 1000x'x' + 1600xNx'.



The configuration of the experimental design as well as the

check point locations are depicted in Figure 4 of Chapter 2

and the observed response values are given in Table 3 of

this chapter.
-1 ^2
A value of the ratio F = [d'V d]/o is calculated at
0 ext
each of the four individual check points (1/3, 1/3, 1/3),

(2/3, 1/6, 1/6), (1/6, 2/3, 1/6), and (1/6, 1/6, 2/3)
"2 ^2
where ext is assumed to have the value et = 144 as

suggested by Kurotori (1966). We also assume without loss

of generality that the degrees of freedom associated
^2
with aet are v = 10. The power of the F test is calculated

at each of the four check points by using the approximation

to the cumulative probabilities of the noncentral F












Table 3. Observed Response Values at the
Pseudocomponent Settings for Kurotori's Rocket Fuel
Experiment--Numerical Example 1.


Observation


Number x'
______ _____


1
2
3
4
5
6
7*
8*
9*
10*


Binder Oxidizer Fuel


Modulus


x 3


0
1/2
1/2
0
1/3
2/3
1/6
1/6


0
1
0
1/2
0
1/2
1/3
1/6
2/3
1/6


0
0
1
0
1/2
1/2
1/3
1/6
1/6
2/3


of Elasticity
Y


2350
2450
2650
2400
2750
2950
3000
2690
2770
2980


*Check Points.


distribution given by Johnson and Kotz (1970, p. 197) to

evaluate


Power = P{FI' ,;,
11, 10;Xi


F.05;1,10}


2 ^2
where 1 = A 8 23/2 ext.
I 1 123 ex


The value of


A = (X* X*A)'V0 (X* X*A) is calculated for each check

point using the values of X*, X*, v and the value of

A = (X'X) X'X2 which is fixed by the {3,2} simplex lattice

design. Since the {3,2} simplex lattice consists of points







only on the boundaries of the triangle (and therefore at

each point at least one of the x! values is equal to zero),
1
then X2 = 0 and A = 0. From the true special cubic model,

123 = 6150.

The calculated value of F as well as the approximate

value for the power at each of the four check points is

given in Table 4. The check point (1/3, 1/3, 1/3) produced

the highest power of the four check points investigated,

supporting the previous results of Example 1 in Section

3.4.1 where (1/3, 1/3, 1/3) was selected as the check point

location with the maximum power when a second order

canonical polynomial was fitted using the {3,2} simplex

lattice design, but the true model was assumed to be special

cubic. Additional support for the point (1/3, 1/3, 1/3)

being optimal is given by the contour plot of values of A1

in Figure 6(d). The highest values of A1 appear near the

centroid (1/3, 1/3, 1/3) where high A1 values translate into
2 2
high X1 values, since Xi = A a123/22 which in turn implies

high power since we know the power is an increasing function

of X1.

As a second part of this example the power of the F

test that is obtained when three replicates are taken at

(1/3, 1/3, 1/3) is compared to the power of the F test that

is obtained when one replicate is taken at the three check

points (2/3, 1/6, 1/6), (1/6, 2/3, 1/6), and (1/6, 1/6,

2/3). These latter three point locations were suggested by

Kurotori for testing lack of fit of his fitted special

















































o 0 0 0
-4


r-1 -4- r- -4
,- -4 ,- O


OC










1-1




mm







M N\N
i- '

- '--


0 *


0w






O


.JU
Q .H

-3


Q)
4J-)
U

-4

0
U





0
m-
4a
4J






(0
>1
cu






4Ja



CO
















00
S0O

-Cu

A 0


u r-
0 0







-c




o




4 r-4








UO-4








A U( -
S*
I -4 4 J
U K M


1-4
o
o
o







o
























r-I- CN


, q --



CN -.4 ,-I
(m)











cMNrN
S S


0 -
4J -4




Z 1







cubic model. The result of this comparison, see Table 4, is

that the three replicates at (1/3, 1/3, 1/3) produce the

test with greater power which again supports the findings of

Example 1 of Section 3.4.1.

All of the check point locations listed in Table 4

produce very high power values (> .999) which is due in part

to the high value of 123 (123 = 6150). If S1 were of
123 123 6123
lower magnitude, then the three replicates at (1/3, 1/3,

1/3) would show a still greater superiority in the power

value compared to the power using the other check points.

This superiority is demonstrated in Table 5 where values of

8123 are listed as 3000 and 1500 and the comparative power

values are listed as 0.998 compared to 0.795 and 0.662

compared to 0.249, respectively. Table 5 also demonstrates

the superior power value for the point (1/3, 1/3, 1/3) when

123 = 3000 or 123 = 1500 and each of the four check points

is used one at a time.

Finally, (1/3, 1/3, 1/3) being the optimal check point

location is seen in Figure 6(c), where contour plots of the

expected difference in the heights of the surfaces are

drawn. The differences in the heights are found by

subtracting the estimated height of the surface obtained

with the fitted second order model from the estimated height

of the surface obtained with the true special cubic model.

The expected difference between the true and fitted surfaces

approaches a maximum the closer one moves to the centroid of

the simplex factor space, so that the optimal check point







89







<4 0O CD CD C
a) m Co C MN C N c
C 0 *

U a4




-I COD0 CO
U- co co 1 o CD m

0
a4


41J -4

C a) m co o o co m a) 4)
0 (a D aN O Oo *10 10 m r- En

C >1* *





0
U a)


0 a U1 m a) 0T co
m 2' 0 0



'0 Z O*o r- *- r*- *< u
Z ; 0 O
Oi 0
*-q 4J V





4 L 0 4 0 *



q m o .m .
0 )








En *' m

a) 0








0 Cr rn)Cr)Cr)

^ N. N. N.. N. 'NN.'N N

a) N N N N NNN N.N U
-; -4 (N i- -4 C '4 M |-4 -4
U .4'













x :I
I


X =1 x= I
2 3
(a) True special cubic surface.








X =I


x :1 x :1
2 3
(c) Expected difference between the
true special cubic surface and
the fitted second order surface.


x2-
(b) Fitted second order surface.


x =1 X :1
2 3

(d) A( X -X XA)' V'( X- X*A)


Figure 6. Contour plots for Numerical Example 1.







location (1/3, 1/3, 1/3) coincides with the point where the

expected difference between the true special cubic surface

and the fitted second order surface is maximum.

Numerical Example 2. In this second numerical example,

we investigate the power of the F test for detecting lack of

fit when a second order canonical polynomial model is fitted

in a mixture system which is in truth represented by a

special cubic model. The true model is assumed to be



E(Y) = 2350x1 + 2450x2 + 2650x3



+ 1000x x3 + 1600x2x3 + 6150x 2X3



which is used to generate hypothetical response observations

at the seven points of the q = 3 simplex centroid design as

well as at three check points. The values of the response

are obtained by adding the value of a pseudorandom normal

variate with mean 0 and variance 144 to each true predicted

response value. The data are given in Table 6.

The response values at the seven points of the simplex

centroid design are used in the least squares normal

equations to obtain the fitted second order model



Y = 2341x1 + 2438x2 + 2630x3


+ 310x1x2 + 1304x x3 + 1970x2x3











Table 6. Generated Response Values--Numerical Example 2.

xl x2 x3 Y

1 0 0 2357
0 1 0 2454
0 0 1 2646
1/2 1/2 0 2403
1/2 0 1/2 2747
0 1/2 1/2 2962
1/3 1/3 1/3 3013
1/3 1/3 1/3 2993
2/3 1/6 1/6 2693
.02 .93 .05 2550

*Check points.









which is to be tested for lack of fit using the test
-1
statistic F = d'V0 d/MSE. The F statistic will be evaluated

at each of the three check points (1/3, 1/3, 1/3),

(2/3, 1/6, 1/6), and (.02, .93, .05), taken one at a time,

and the power of the test at the three check point locations

will be calculated and compared. The test is lower tailed

for all check point locations (since R = A1 A2 is negative

for all check point locations) and thus the power is defined

as



PI" ( F;, }.
1,1;XIX2 .95;1,1




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs