5.
Process Improvement
5.5. Advanced topics 5.5.3. How do you optimize a process? 5.5.3.1. Single response case


Regions where quadratic models or even cubic models are needed occur in many instances in industry 
After a few steepest ascent (or descent) searches, a firstorder model
will eventually lead to no further improvement or it will exhibit lack
of fit. The latter case typically occurs when operating conditions
have been changed to a region where there are quadratic (secondorder)
effects present in the response. A secondorder polynomial can be
used as a local approximation of the response in a small region where,
hopefully, optimal operating conditions exist. However, while a
quadratic fit is appropriate in most of the cases in industry, there
will be a few times when a quadratic fit will not be sufficiently
flexible to explain a given response. In such cases, the analyst
generally does one of the following:


Procedure: obtaining the estimated optimal operating conditions  
Second order polynomial model 
Once a linear model exhibits lack of fit or when significant curvature
is detected, the experimental design used in Phase I (recall that a
2^{kp} factorial experiment might be
used) should be augmented with axial runs on each factor to form what
is called a central composite design. This experimental design
allows estimation of a secondorder polynomial of the form


Steps to find optimal operating conditions 
If the corresponding analysis of variance table indicates no lack of
fit for this model, the engineer can proceed to determine the estimated
optimal operating conditions.


Chemical experiment example  We illustrate these steps using the chemical experiment discussed previously. For a technical description of a formula that provides the coordinates of the stationary point of the surface, see Technical Appendix 5C.  
Example: Second Phase Optimization of Chemical Process  
Experimental results for axial runs 
Recall that in the chemical experiment, the
ANOVA table,
obtained from using an experiment run around the coordinates
X_{1} = 189.5, X_{2} = 350,
indicated significant curvature effects. Augmenting the
2^{2} factorial experiment with axial runs at
\( \pm \alpha = \pm \sqrt{2} \)
to achieve a rotatable central composite experimental design, the
following experimental results were obtained:


ANOVA table 
The ANOVA table corresponding to a cubic model with an interaction term (contained
in the quadratic sumofsquares partition) is
SUM OF MEAN F SOURCE SQUARES DF SQUARE VALUE PROB > F MEAN 51418.2 1 51418.2 Linear 1113.7 2 556.8 5.56 0.024 Quadratic 768.1 3 256.0 7.69 0.013 Cubic 9.9 2 5.0 0.11 0.897 RESIDUAL 223.1 5 44.6 TOTAL 53533.0 13 

Lackoffit tests and auxillary diagnostic statistics 
From the ANOVA table, the linear and quadratic effects are significant. The
lackoffit tests and auxiliary diagnostic statistics for linear, quadratic,
and cubic models are:
SUM OF MEAN F MODEL SQUARES DF SQUARE VALUE PROB > F Linear 827.9 6 138.0 3.19 0.141 Quadratic 59.9 3 20.0 0.46 0.725 Cubic 49.9 1 49.9 1.15 0.343 PURE ERROR 173.2 4 43.3 ROOT ADJ PRED MODEL MSE RSQR RSQR RSQR PRESS Linear 10.01 0.5266 0.4319 0.2425 1602.02 Quadratic 5.77 0.8898 0.8111 0.6708 696.25 Cubic 6.68 0.8945 0.7468 0.6393 3466.71The quadratic model has a larger pvalue for the lack of fit test, higher adjusted R^{2}, and a lower PRESS statistic; thus it should provide a reliable model. The fitted quadratic equation, in coded units, is


Step 1:  
Contour plot of the fitted response function 
A contour plot of this function (Figure 5.5) shows that it appears to
have a single optimum point in the region of the experiment (this
optimum is calculated below to be (0.9285, 0.3472), in coded
x_{1}, x_{2} units, with a
predicted response value of 77.59).


3D plot of the fitted response function 
Since there are only two factors in this example, we can also obtain a
3D plot of the fitted response against the two factors (Figure 5.6).
FIGURE 5.6: 3D Plot of the Fitted Response in the Example 

Step 2:  
Optimization point  An optimization routine was used to maximize \( \hat{Y} \). The results are \( X_{1}^{*} = 161.64^{\circ}C \), \( X_{2}^{*} = 367.32 \) minutes. The estimated yield at the optimal point is \( \hat{Y}(X^{*}) = 77.59 \)%,  
Step 3:  
Confirmation experiment  A confirmation experiment was conducted by the process engineer at settings X_{1} = 161.64, X_{2} = 367.32. The observed response was \( \hat{Y}(X^{*}) = 76.5 \)%, which is satisfactorily close to the estimated optimum.  


Technical Appendix 5C: Finding the Factor Settings for the Stationary Point of a Quadratic Response  
How to find the maximum or minimum point for a quadratic response 


Nature of the stationary point is determined by B  The nature of the stationary point (whether it is a point of maximum response, minimum response, or a saddle point) is determined by the matrix B. The twofactor interactions do not, in general, let us "see" what type of point x^{*} is. One thing that can be said is that if the diagonal elements of B (b_{ii}) have mixed signs, x^{*} is a saddle point. Otherwise, it is necessary to look at the characteristic roots or eigenvalues of B to see whether B is "positive definite" (so x^{*} is a point of minimum response) or "negative definite" (the case in which x^{*} is a point of maximum response). This task is easier if the twofactor interactions are "eliminated" from the fitted equation as is described in Technical Appendix 5D.  
Example: computing the stationary point, Chemical Process experiment  
Example of computing the stationary point 
The fitted quadratic equation in the chemical experiment discussed in
Section 5.5.3.1.1 is, in coded units,
and


Technical Appendix 5D: "Canonical Analysis" of Quadratic Responses  
Case for a single controllable response  Whether the stationary point X^{*} represents a point of maximum or minimum response, or is just a saddle point, is determined by the matrix of secondorder coefficients, B. In the simpler case of just a single controllable factor (k=1), B is a scalar proportional to the second derivative of \( \hat{Y}(x) \) with respect to x. If \( d^{2}\hat{Y}/dx^{2} \) is positive, recall from calculus that the function \( \hat{Y}(x) \) is convex ("bowl shaped") and x^{*} is a point of minimum response.  
Case for multiple controllable responses not so easy  Unfortunately, the multiple factor case (k>1) is not so easy since the twofactor interactions (the offdiagonal elements of B) obscure the picture of what is going on. A recommended procedure for analyzing whether B is "positive definite" (we have a minimum) or "negative definite" (we have a maximum) is to rotate the axes x_{1}, x_{2}, ..., x_{k} so that the twofactor interactions disappear. It is also customary (Box and Draper, 1987; Khuri and Cornell, 1987; Myers and Montgomery, 1995) to translate the origin of coordinates to the stationary point so that the intercept term is eliminated from the equation of \( \hat{Y}(x) \). This procedure is called the canonical analysis of \( \hat{Y}(x) \).  
Procedure: Canonical Analysis  
Steps for performing the canonical analysis 


Eigenvalues that are approximately zero 
If some λ_{i}
≈ 0, the fitted ellipsoid


Canonical analysis typically performed by software  Software is available to compute the eigenvalues λ_{i} and the orthonormal eigenvectors e_{i}; thus there is no need to do a canonical analysis by hand.  
Example: Canonical Analysis of Yield Response in Chemical Experiment  
B matrix for this example 
Let us return to the chemical experiment
example to illustrate the method. Keep
in mind that when the number of factors is small (e.g., k=2
as in this example) canonical analysis is not recommended in practice
since simple contour plotting will provide sufficient information. The
fitted equation of the model yields


Compute the eigenvalues and find the orthonormal eigenvectors  To compute the eigenvalues λ_{i}, we have to find all roots of the expression that results from equating the determinant of B  λ_{i}I to zero. Since B is symmetric and has real coefficients, there will be k real roots λ_{i}, i = 1, 2, ..., k. To find the orthonormal eigenvectors, solve the simultaneous equations \( (B  \lambda_{i}I)e_{i} = 0 \) and \( e_{i}^{'}e_{i} = 1 \).  
Canonical analysis results 
The results of the canonical analysis are as follows:
Eigenvectors Eigenvalues X1 X2 4.973187 0.728460 0.685089 9.827317 0.685089 0.728460Notice that the eigenvalues are the two roots of
