The purpose of this article is to provide an overview of response surface methodology (RSM) which includes the modeling of a response function, the corresponding choice of design, and the determination of optimum conditions. In addition, the use of RSM in agricultural and food sciences is highlighted by citing several examples taken from a variety of applied journals.

Response surface methodology (RSM) consists of a group of mathematical and statistical techniques concerned with

1. The selection and construction of an appropriate design that can provide adequate and reliable information concerning a certain Response variable, denoted by $y$ .
2. The determination of a suitable model that best fits the data that can be generated from using the design chosen in (1). Such a model gives an approximate functional relationship between the response variable $y$ and a set of Control variables believed (by the experimenter) to have an effect on the response $y$ . These variables are denoted by $x_1,x_2,\mathrm{...},x_k$ .
3. The finding of optimal settings on the control variables that produce maximum (or minimum) response values within a certain region of interest.

From the historical point of view, the article by Box and Wilson¹ is considered to be the first for having laid the foundations for RSM. Some of the early papers that also contributed to the development of RSM include those by Box and Hunter² and Box and Draper.³ Several review papers on RSM were subsequently published starting with the article by Hill and Hunter⁴ that emphasized practical applications of RSM in the chemical and processing fields. This was followed by more recent reviews by Myers et al.^5–7 A review of RSM from a biometric viewpoint was given by⁸ who emphasized biological applications rather than applications in the physical and engineering sciences. In addition, several books were written on the subject of RSM by Myers, Khuri and Cornell, Box and Draper and Myers et al.^9–12

In the early development of RSM, from 1951 up through the 1970s, the main focus of attention was on controlled experiments that can, for the most part, be performed in a laboratory. This was particularly suited in an industrial setting with possible applications in the physical and engineering sciences. The review paper by Hill and Hunter⁴ made reference to several applications in the chemical industry. It should be pointed out here that the seminal work of Box and Wilson¹ occurred at a major chemical company, and Box himself was initially trained as a chemist in England.

One of the basic characteristics of a response surface investigation is its sequential nature whereby experiments are performed in stages. Information acquired from one set of experiments is used to plan the strategy for a follow-up set of experiments. This sequential pattern of experimentation was suggested by Box and Youle.¹³ Such an approach works well in an industrial setting since the response values in a given stage can be obtained in a relatively short time. However, this approach may not be feasible in an agricultural setting where experiments typically require long periods between stages. Furthermore, it is quite common in an agricultural experiment that the results of a single experiment may be the only ones available rather than having a series of experiments obtained sequentially. Another main difference between an industrial experiment and an agricultural one has to do with factors' levels. In an industrial experiment, the levels of quantitative factors can be accurately controlled and measured. This, however, may be difficult to do in an agricultural experiment. Edmondson¹⁴ Outlined several other major distinctions between industrial and agricultural experiments.

Mead and Pike⁸ were the first to bring attention to the role of RSM in agriculture and biometric research in general. They also provided a survey of bio- logical and agricultural journals that used response surface methods. However, they emphasized the use of nonlinear models to describe the behavior of biological data. Such models are usually referred to as Mechanistic. The traditional modeling scheme used in the early work on RSM was based on the so-called Empirical modeling where low-degree polynomials are used to fit the response data. These polynomials are chosen and tested on the basis of the observed data, but are not selected on the basis of information that pertain to certain chemical of physical laws, as is the case with mechanistic models.

While certain approaches used in RSM may not be quite suitable in an agricultural setting, there is a lot to be gained from using certain response surface techniques in an agricultural experiment. The purpose of this chapter is to provide an expose of some basic methods and models used in RSM. Several applications of RSM in agricultural and food sciences will also be mentioned using examples taken from the corresponding literature.

Let $y$ be a response variable of interest, and let $x_1,x_2,\mathrm{...},x_k$ denote control variables believed to have an effect on $y$ . As was mentioned earlier, one of the objectives of RSM is to establish a functional relationship between $y$  and its associated control variables. Such a relationship is, in general, unknown, but can be approximated by a low-degree polynomial model of the form $y=f^{\text{'}}\left(x\right)\beta +\in$ , where $x={\left(x_1,x_2,\mathrm{...},x_k\right)}^{\text{'}}$ $f\left(x\right)$  is a vector function of p elements that consists of powers and cross products of powers of $x_1,x_2,\mathrm{...},x_k$  up to a certain degree denoted $d\left(\ge 1\right),$ $\beta$  is a vector of p unknown constant coefficients referred to as parameters, and _ is a random experimental error assumed to have a zero mean. If model (2.1) provides an adequate representation of the response, then the quantity $f^{\text{'}}\left(x\right)\beta$ represents the so-called Mean response, or the Expected Value of $y$ , and is denoted by $\mu \left(x\right)$ .

Two models are commonly used in RSM; they are the first-degree model (d = 1),

$y={\beta }_0+{\displaystyle \sum _{i=1}^k{\beta }_i}{x}_i+\in$ ;              (2.2)

and the second-degree model (d = 2),

$y={\beta }_0+{\displaystyle \sum _{i=1}^k{\beta }_i}{x}_i+{\displaystyle \sum _{i<j}{\displaystyle \sum {\beta }_{ij}{x}_i{x}_j+{\displaystyle \sum _{i=1}^k{\beta }_{ii}}{x}_i^2+\in }}$ ;         (2.3)

Models (2.2) and (2.3) are special cases of model (2.1). The first model is usually used in the initial phase of the experiment representing part of an exploratory process to assess the important factors. The second model is subsequently developed after a series of experiments have been performed resulting in the identification of the important factors to be considered in the experiment. At this stage, the experimenter will be ready to use the model in doing data analysis to determine significance of the model's parameters, estimate the mean response, and to arrive at optimum operating conditions on the control variables that result in a maximum or a minimum response over a certain region of interest which we denote by $R$ .

In order to estimate the unknown parameters in model (2.1), a series on $n\left(>p\right)$  experiments are performed in each of which the response y is measured for specified settings of the control variables. The totality of these settings constitutes the so called response surface design, or just design, which can be represented by a matrix, denoted by $D$ , of order $n\times k$ called the design matrix,

$D=\left[\begin{array}{cccc}{x}_{11}& x_{12}& \mathrm{...}& x_{1k}\\ x_{21}& x_{22}& \mathrm{...}& x_{2k}\\ .& .& \mathrm{...}& .\\ .& .& \mathrm{...}& .\\ .& .& \mathrm{...}& .\\ x_{n1}& x_{n2}& \mathrm{...}& x_{nk}\end{array}\right]$ ;(2.4)

Where $x_{ui}$ denotes the $u^{th}$  design setting of $x_i\left(i=1,2,\mathrm{...},k;u=1,2,\mathrm{...},n\right)$ . Each row of $D$  represents a point, referred to as a design point, in a k-dimensional ${\in }_u$  clidean space. The response value obtained at the $u^{th}$  setting of $x$ , namely   $x_u={\left(x_{u1},x_{u2},\mathrm{...},x_{uk}\right)}^{\text{'}}$ $\left(u=1,2,\mathrm{...},n\right)$ , is denoted by $y_u$ . From (2.1) we thus have

$y_u=f^{\text{'}}\left(x_u\right)\beta +{\in }_u,$ ; $u=1,2,\mathrm{....},n$ (2.5)

Where ${\in }_u$ denotes the error term at the $u^{th}$ experimental run. Model (2.5) can be expressed in matrix form as

$y=X\beta +\in$ ;(2.6)

Where $y={\left(y_1,y_2,\mathrm{...},y_n\right)}^{\text{'}}$ 0, $X$  is a matrix of order $n\times p$ whose $u^{th}$  row is $f^{\text{'}}\left(x_u\right)$ , and $\in ={\left({\in }_1,{\in }_2,\mathrm{....},{\in }_n\right)}^{\text{'}}$ . Note that the first column of $X$  is the column of ones, $1_n$ .

A common method for estimating the parameter vector in (2.6) is the one based on the method of ordinary least squares. This method requires that has a zero mean and a variance-covariance matrix given by ${\sigma }^2{I}_n$  [see, for example, Khuri and Cornell (Section 2.3)].¹⁰ In this case, the least-squares estimator of $\beta$ , denoted by $\widehat{\beta }$ , is given by

$\widehat{\beta }={\left(X^{\text{'}}X\right)}^{-1}{X}^{\text{'}}y$ ; (2.7)

The variance-covariance matrix of $\widehat{\beta }$ is then of the form

$Var\left(\widehat{\beta }\right)={\left(X^{\text{'}}X\right)}^{-1}{X}^{\text{'}}\left({\sigma }^2{I}_n\right)X{\left(X^{\text{'}}X\right)}^{-1}$

$={\sigma }^2{\left(X^{\text{'}}X\right)}^{-1}$ :        (2.8)

An estimate $\widehat{\mu }\left(x_u\right)$ , of the mean response at $x_u$ , is obtained by replacing $\beta$  by $\widehat{\beta }$ , that is,

  $\widehat{\mu }\left(x_u\right)=f^{\text{'}}\left(x_u\right)\widehat{\beta },$ $$u=1,2,\mathrm{...},n.$$ :       (2.9)

The quantity $f^{\text{'}}\left(x_u\right)\widehat{\beta }$ also gives the value of the predicted response, $\widehat{y}\left(x_u\right)$ , at the $u^{th}$  design point $\left(u=1,2,\mathrm{...},n\right)$ . In general, at any point, $x$ , in the region $R$ , the predicted response $\widehat{y}\left(x\right)$  is

$\widehat{y}\left(x\right)=f^{\text{'}}\left(x\right)\widehat{\beta },$ $x\in R$ : (2.10)

Using formula (2.8), the variance of $\widehat{y}\left(x\right)$ is then given by [Khuri and Cornell (1996, Section 6.1) or Khuri (2009, Section 12.4)]

$Var\left[\widehat{y}\left(x\right)\right]={\sigma }^2{f}^{\text{'}}\left(x\right){\left(X^{\text{'}}X\right)}^{-1}f\left(x\right)$ : (2.11)

The expression on the right-hand side of (2.11) is called the Prediction variance. This is an important quantity since the quality of prediction depends on the size of this variance. Also, the determination of optimum operating conditions on the control Variables requires that the prediction variance be small over the region of interest, $R$ , in order to arrive at reliable information about the true optimum of the response over $R$ . This is of course dependent on the assumption that the postulated model in (2.1) does not suffer from Lack of fit (LOF). For a study of LOF of a fitted response surface model, see, for example, Khuri and Cornell (Section 2.6).¹⁰

The choice of the design matrix, $D$ , is quite important in any response surface investigation since the prediction variance depends on $D$  as can be seen from formula (2.11). Some common design properties are

1. Orthogonality

A design $D$  is said to be orthogonal if the matrix $X^{\text{'}}X$  is diagonal, where $X$ is the model matrix in (2.6). In this case, the elements of $\widehat{\beta }$  will be uncorrelated since the off-diagonal elements of $Var\left(\widehat{\beta }\right)$ in (2.8) will be zero. If the error vector $\in$ in (2.6) is assumed to be normally distributed as $N\left(0,{\sigma }^2{I}_n\right)$ , then these elements will be also stochastically independent. This makes it easier to test the significance of the unknown parameters in the model.

2. Rotatability

A design $D$ is said to be rotatable if the prediction variance in (2.11) is constant at all points that are equidistant from the design center, which, by a proper coding of the control variables, can be chosen to be the point at the origin of the $k$ -dimensional coordinates system. It follows that $Var\left[\widehat{y}\left(x\right)\right]$ is constant at all points that fall on the surface of a hyper sphere centered at the origin, if the design is rotatable. This causes the prediction variance to remain unchanged under any rotation of the coordinate axes. In addition, if optimization of $\widehat{y}\left(x\right)$ is desired on concentric hyper spheres, as in the application of ridge analysis, which will be discussed later, then it would be desirable for the design to be rotatable. This makes it easier to compare the values of $\widehat{y}\left(x\right)$ on a given hyper sphere as all such values will have the same variance.

The necessary and sufficient condition for a design to be rotatable was given by Box and Hunter.² See also Appendix 2B in^15,16 introduced a measure that quantities the amount of rotatability in a response surface design. This measure can be helpful in comparing designs on the basis of rotatability, assessing the extent of departure from rotatability, and in improving readability by a proper augmentation of a non rotatable design.

3. Optimal designs

Optimal designs are those that are constructed on the basis of a certain optimality criterion that pertains to the 'closeness' of the predicted response, $\widehat{y}\left(x\right)$ , to the mean response, $\mu \left(x\right)$ over a certain region of interest, $R$ . The design criteria that address the minimization of the variance associated with the estimation of model (2.1)'s unknown parameters are called Variance-related criteria. The most prominent of such criteria is the D-Optimality criterion which maximizes the determinant of the matrix $X^{\text{'}}X$ . This amounts to the minimization of the size of the confidence region on the vector $\beta$ in model (2.6). Another variance-related criterion that is closely related to D-optimality is the G-Optimality criterion which requires the minimization of the maximum over $R$  of the prediction variance in (2.11).

Such variance-related criteria are often referred to as Alphabetic optimality. They are meaningful when the fitted model in (2.1) represents the true relation- ship connecting $y$  to its control variables. There are, however, many situations Where this is not the case due to fitting the "wrong model". This results in the So-called Model bias. For example, a first-degree model may be fitted to a data set when in reality a second-degree model would be a better representative of the response data. Box and Draper³ placed a great emphasis on the role of bias in the choice of a response surface design and advocated choosing designs that minimized model bias. They considered the minimization of bias to be a very important design criterion, and in certain cases, even more important than the variance-related criteria.

Designs for first-degree models

Designs for fitting first-degree models as in (2.2) are called first-order designs. The most common of such designs are the $2^k$ factorial ( $k$  is the number of control variables), and the Plackett-Burman design.

The 2k Factor $2^k$ ial Design: In a  factorial design, each control variable is measured at two levels, which can be coded to take the values, -1, 1, that correspond to the so-called low and high levels, respectively, of each variable. This design consists of all possible combinations of such levels of the $k$ factors. Thus, each row of the design matrix $D$ in (2.4) consists of all 1's, all -1's, or a combination thereof. It therefore represents a particular treatment combination. In this case, the number, $n$ , of experimental runs is equal to $2^k$  provided that no single design point is replicated more than once. For example, in an agricultural experiment, two levels of fertilizer A are combined with two levels of fertilizer B in order to study their effects on the yield of a certain vegetable crop over a certain period of time. This results in a $2^2$ factorial experiment with four treatment combinations.

If $k$ is large $\left(k\ge 5\right)$ , the $2^k$  design requires a large number of design points. In this case, fractions of $2^k$  can be considered. For example, we can consider a one-half- fraction design which consists of one-half the number of points of a $2^k$  design, or a one-fourth-fraction design which consists of one- fourth the number of points of a $2^k$  design. In general, a $2^{-m}th$  fraction of a $2^k$  design consists of $2^{k-m}$  points from a full $2^k$  design. Here, m is a positive integer such that $2^{k-m}\ge k+1$ so that all the $k+1$  parameters in model (2.2) can be estimated. The construction of fractions of a $2^k$  design is carried out in a particular manner, a description of which can be found in several experimental design textbooks, such as Box & Raktoe, et al.^17–19 See also Chapter 3 in Khuri and Cornell.¹⁵

The plackett-burman design: The Plackett-Burman design allows two levels for each of the $k$  control variables, just like a $2^k$ design, but requires a much smaller number of experimental runs, especially if $k$  is large. It is therefore more economical than the $2^k$  design. Its number, $n$ , of design points is equal to $k+1$ , which is the same as the number of parameters in model (2.2). In this respect, the design is said to be saturated because its number of design points is equal to the number of parameters to be estimated in the model. Furthermore, this design is available only when n is a multiple of 4. Therefore, it can be used when the number, $k$ , of control variables is equal to 3, 7, 11, 15, ....

To construct a Plackett-Burman design in $k$  variables, a row is first selected whose elements are equal to -1 or 1 such that the number of 1's is $\frac{k+1}{2}$ and the number of -1's is $\frac{k-1}{2}$ . The next $k-1$ rows are generated from the first row by shifting it cyclically one place to the right $k-1$ times. Then, a row of negative ones is added at the bottom of the design. For example, for $k=7$ , the design matrix, $D$ , has 8 points whose coordinates are $x_1,x_2,\mathrm{...},x_{7,}$ , and is of the form

$D=\left[\begin{array}{ccccccc}1& 1& 1& -1& 1& -1& -1\\ -1& 1& 1& 1& -1& 1& -1\\ -1& -1& 1& 1& 1& -1& 1\\ 1& -1& -1& 1& 1& 1& -1\\ -1& 1& -1& -1& 1& 1& 1\\ 1& -1& 1& -1& -1& 1& 1\\ 1& 1& -1& 1& -1& -1& 1\\ -1& -1& -1& -1& -1& -1& -1\end{array}\right]$ :

Design arrangements for $k=3,7,11,\mathrm{...},99$ factors can be found in Plackett –Burman (1946).

Designs for second-degree models

These are designs for fitting second-degree models as in (2.3), which has $p=1+2k+\frac{1}{2}k\left(k-1\right)$ parameters (they are also referred to as second-order designs). The number of distinct design points of such design must therefore be at least equal to p. The design settings are sometimes coded so that $\frac{1}{n}{\displaystyle {\sum }_{u=1}^n{x}_{ui}=0}$ and
$\frac{1}{n}{\displaystyle {\sum }_{u=1}^n{x}_{ui}^2=1},$ , $i=1,2,\mathrm{...},k,$ , where $n$  is the number of experimental runs and $x_{ui}$  is the $u^{th}$  setting of the $i^{th}$ control variable $\left(u=1,2,\mathrm{...},n,i=1,2,\mathrm{...},k\right)$ .

The most frequently-used second-order designs are the $3^k$  factorial, central com- posite, and the Box-Behnken designs.

The   $3^k$ Factorial Design: The  factorial design consists of all the combinations of the levels of the $k$ control variables which have three levels each. If the levels are equally spaced, then they can be coded so that they correspond to -1, 0, 1. For example, for $k=2$ , the design matrix, in coded form, consists of 9 points as shown below:

$D=\left[\begin{array}{cc}-1& -1\\ -1& 0\\ -1& 1\\ 0& -1\\ 0& 0\\ 0& 1\\ 1& -1\\ 1& 0\\ 1& 1\end{array}\right]$ :

The number of design points for a $3^k$  design is $3^k$ , which can be very large for a large $k$ . The use of full factorial designs is therefore not feasible if the number of experimental units is limited. Fractions of a $3^k$  design can be considered to reduce the cost of running such an experiment. A general procedure for constructing fractions of $3^k$ is described, for example, in.¹⁹ See also McLean and Anderson.²⁰

The central composite design: The central composite design (CCD) is perhaps the most popular of all second-order designs. It was first introduced in¹ as an alternative to the design. This design consists of the following three parts:

1.  A complete (or a fraction of) $2^k$  factorial design whose factors' levels are coded as -1, 1. This is called the factorial portion of the design.
2.  An axial portion consisting of $2k$  points arranged so that two points are chosen on the axis of each control variable at a distance of $\alpha$  from the design center (chosen as the point at the origin of the coordinates system). We refer to $\alpha$ as the axial parameter.
3.  A certain number, $n_0$ , of replications at the design center $\left(n_0\ge 1\right)$ . This is called the center-point portion.

 Thus, the total number of design points in a CCD is $n=2^k+2k+n_0$ . For example, a CCD for $k=2,\alpha =\sqrt{2},$ $n_0=2$ has the form

$D=\left[\begin{array}{cc}-1& -1\\ 1& -1\\ -1& 1\\ 1& 1\\ -\sqrt{2}& 0\\ \sqrt{2}& 0\\ 0& -\sqrt{2}\\ 0& \sqrt{2}\\ 0& 0\\ 0& 0\end{array}\right]$ :

We note that the CCD is obtained by augmenting a first-order design, namely, the and then $n_0$  center-point replications. This design is usually developed in a manner consistent with the sequential nature of a response surface investigation in starting with a first- $2^k$  factorial with additional experimental runs, namely, the $2k$  axial points order design, to fit a first-degree model, followed by the addition of design points to _t the larger second-degree model. The first-order design serves in a preliminary phase to get initial information about the response system and to assess the importance of the factors in a given experiment. The additional experimental runs are chosen for the purpose of estimating all the $p$ parameters in model (2.3). The fitted model is then used in the determination of optimum operating conditions on the control variables over the region of experimentation.

When $k$  is large $\left(k\ge 5\right)$ , the factorial portion can be replaced by a fraction of a $2^k$  design. For example, for $k=5$ , a one-half fraction of $2^5$  can be used giving a total of 16 points in the factorial portion instead of 32 (for more details about fractionating in the factorial portion, see Khuri and Cornell, 1996, Section 4.5.3).

A CCD can become rotatable by assigning the value $F^{1/4}$ to the axial parameter, $\alpha$ , where $F$ denotes the number of points in its factorial portion, that is, $\alpha =F^{1/4}$ . In addition, the number of center-point replications, $n_0$ , can be chosen so that a rotatable CCD will have the additional property of orthogonality (see Section 3). Note that orthogonality of a second-order design is attainable only after expressing model (2.3) in terms of orthogonal polynomials as explained in.² See also.¹⁵ In particular, Table 4.3 in Khuri and Cornell's book can be used to determine the value of $n_0$ in order for a rotatable CCD to have the additional orthogonality property.

The box-behnken design: This design was developed by Box GEP et al.²¹ It provides three levels for each factor and consists of a particular subset of the factorial combinations from the factorial design. The actual construction of such a design is described in the three RSM books,¹¹ Khuri and Cornell (1996, Section 4.5.2), and.¹² Some Box-Behnken designs are rotatable, but this design is not always rotatable. Box GEP²¹ list a number of design arrangements for $k=3,4,5,6,7,9,10,11,12,and16$  control variables.

The san cristobal design: Rojas²² introduced a variant of the CCD, called the San Cristobal Design (SCD), for sugar farming experiments. It is utilized in situations where the levels of $k$ control variables are restricted to be nonnegative, as is the case with fertilizers experiments. The SCD consists of $2^k$ factorial points combined with center and axial points, all contained within the positive orthant. It also includes a control where no fertilizers are applied. More recently, the performance of this design was evaluated by Haines LM²³ who reviewed some of its properties.

None.

None.

1. Box GEP, Wilson KB. On the experimental attainment of optimum conditions. Journal of the Royal Statistical Society. 1951;B 13:1–45.
2. Box GEP, Hunter JS. Multifactor experimental designs for ex-ploring response surfaces. Ann Math Statist. 1957;28(1):195–241.
3. Box GEP, Draper NR. A basis for the selection of a response surface design. J Amer Stat Assoc. 1959;54:622–654.
4. Hill WJ, Hunter WG. A review of response surface methodology: A literature survey. Technometrics. 1966;8(4):571–590.
5. Myers RH, Khuri AI, Carter WH. Response surface methodology: 1966-1988. 1989;31:137–157.
6. Myers RH. Response surface methodology - Current status and future directions. J Qual Technology. 1999;31:30–44.
7. Myers RH, Montgomery DC, Vining GG, et al. Response surface methodology: A retrospective and literature survey. J Qual Technology. 2004;36:53–77.
8. Mead, R, Pike DJ.A review of response surface methodology from a biometric viewpoint. Biometrics. 1975;31:803–851.
9. Myers RH. Response Surface Methodology. Allyn and Bacon, Boston; 1971.
10. Khuri, AI. Multiresponse surface methodology. In: Handbook of Statistics, Vol. 13, Subir Ghosh, CR Rao, editors. Elsevier Science B. V, Amsterdam; 1996. p. 377–406.
11. Box GEP, Draper NR. Response Surfaces, Mixtures, and Ridge Analyses, 2nd ed. 2007.
12. Myers RH, Montgomery DC, Anderson-Cook CM. Response Surface Methodology. Wiley, Hoboken, New Jersey; 2009.
13. Box GEP, Youle PV. The exploration and exploitation of response surfaces: An example in the link between the fitted surface and the basic mechanism of the system Biometrics.1955;11(3):287–323.
14. Edmondson RN. Agricultural response surface experiments based on four-level factorial designs. Biometrics. 1991;47(4):1435–1448.
15. Khuri AI, Cornell JA. Response Surfaces, 2nd ed. Dekker, New York, USA; 1996.
16. Khuri AI. A measure of rotatability for response surface designs, Techno-metrics.1988;30:95–104.
17. Box GEP, Hunter WG, Hunter JS. Statistics for Experimenters; 1978.
18. Ratkoe BL, Hedayat A, Federer WT. Factorial Designs. Wiley, New York, USA; 1981.
19. Montgomery DC. Design and Analysis of Experiments. 6th ed. Wiley, New York, USA; 2005.
20. McLean RA, Anderson VL. Applied Factorial and Fractional De-signs. Dekker, New York; 1984.
21. Box GEP, Behnken DW. Some new three-level designs for the study of quantitative variables. Technometrics. 1960;2(4):455–475.
22. Rojas BA. The San Cristobal design for fertilizer experiments, Proceedings of the 11th Congress of International Society for Sugar Cane Technology. 1962. p. 197–203.
23. Haines LM. Evaluating the performance of non-standard designs: The San Cristobal design. In: Response Surface Methodology and Related Topics. In: Andre I Khuri, editor. World Scientific, Singapore; 2006. p. 251–281.
24. Hoerl AE. Optimum solutions of many variable equations. Chem Eng Prog. 1953;55:69–78.
25. Draper NR. Ridge analysis of response surfaces. Technometrics. 1963;5(4):469–479.
26. Khuri AI and Myers RH. Modified ridge analysis. Technometrics. 1979;21:467–473.
27. SAS Institute, Inc. SAS Online Doc, Version 9.1.3. Author, Cary, North Carolina; 2008.
28. Tseo, CL, Deng JC, Cornell JA, et al. Effect of washing treatment on quality of minced mullet flesh. J Food Sci. 1983;48:163–167.
29. Schmidt RH, Illingworth BL, Deng JC, et al. Multiple regression and response surface analysis of the effects of calcium chloride and cysteine on heat-induced whey protein gelation. J Agric Food Chem. 1979;27:529–532.
30. Khuri AI, Conlon M. Simultaneous optimization for multiple responses represented by polynomial regression functions. Technometrics. 1981;23:363–375.
31. Khuri AI. The analysis of multiresponse experiments: A review. In: Statistical Design and Analysis of Industrial Experiments. Subir Ghosh, editor. Dekker, New York, USA; 1990. p. 231–246.
32. Evans RA, Easi DA, Book DN, et al. Quadratic response surface analysis of seed-germination trials. Weed Science. 1982;30(4):411–416.
33. Keisling TC, Oliver LR, Crowley RH, et al. Potential use of response surface analyses for weed management in soybeans (Glycine max), Weed Sci. 1984;32:552–557.
34. Broudiscou LP, Papon Y, Broudiscou AF. Effects of minerals on feed degradation and protein synthesis by rumen micro-organisms in a dual e uent fermenter. Reprod Nutr Dev. 1999;39:255–268.
35. Salawu IS, Adeyemi RA, Aremu TA. Modified inverse polynomial and ordinary polynomial as a response surface model: A case study of nitrogen, phosphate and potassium levels on the yield of maize. Inter J Pure Appld Sciences. 2007;1:18–24.
36. Nelder JA. Inverse polynomials, a useful group of multi-factor response function. Biometrics.1962;22:128–141.
37. Diniz FM, Martin AM. Use of response surface methodology to describe the combined effects of pH, temperature, and E/S ratio on the hydrolysis of dog sh (Squalus acanthias) muscle. Inter J Food Sci Tech. 2003;31:419–426.
38. Cao W, Zhang C, Hong P, et al. Optimization of the production of shrimp waste protein hydrolysate using microbial proteases adopting response surface methodology. Inter J Food Sci Tech. 2009;44:1602–1608.
39. Jiang G, Jiang Y, Yang B, et al. Structural characteristics and antioxidant activities of oligosaccharides from longan fruit pericarp. J Aqric Food Chem. 2009;57(19):9293–9298.
40. Zhang DH, Bai S, Dong XY, Sun Y. Optimization of lipase-catalyzed regioselective acylation of pyridoxine (vitamin B6). J Agric Food Chem. 2007;55:4526–4531.
41. Bas D, Boyaci IH. Modeling and optimization I: Usability of response surface methodology. Journal of Food Engineering. 2007;78(3):836–845.
42. Box GEP, Draper NR. The choice of a second-order rotatable design. Biometrika.1965;50(1–2):335–352.
43. Khuri, AI. Linear Model Methodology. Chapman & Hall/CRC, Boca Raton, Florida; 2009.
44. Plackett RL, Burman JP. The design of optimum multifactorial experiments. Biometrika. 1946;33:305–325.