Question

1. A chemical engineer is interested in determining the operating conditions that maximize the yield of a process. The data is with coded variables (x1, x2)
2. a). Read in the dataset into R (Write your R code).

           b). Apply Response Surface Method in R, use first-order model (Write your R

            code).

           c). Is the first-order model fit the data? (Write down the P-value for Lack of fit

                test, and answer the question at _ = 0.05).

          d). Apply Response Surface Method in R, use second-order model (Write your R

                code). Write down the stationary point (if exists).

x1	x2	Response
-1	-1	76.5
-1	1	77
1	-1	78.1
1	1	79.5
0	0	79.9
0	0	80.3
0	0	80
0	0	79.7
0	0	79.8
1.414	0	78.4
-1.414	0	75.6
0	1.414	78.5
0	-1.414	77.2

Answer 1

Hi

I think this will help you and make you understand the logic of it.

If not please feel free to comment.

I am happy to take your comment.

> library(rsm)

> x1=c(-1,-1,1,1,0,0,0,0,0,1.414,-1.414,0,0)

> x2=c(-1,1,-1,1,0,0,0,0,0,0,0,1.414,-1.414)

> response=c(76.5,77,78.1,79.5,79.9,80.3,80,79.7,79.8,78.4,75.6,78.5,77.2)

> d1=data.frame(x1,x2,response)

> d1

       x1     x2 response

1  -1.000 -1.000     76.5

2  -1.000  1.000     77.0

3   1.000 -1.000     78.1

4   1.000  1.000     79.5

5   0.000  0.000     79.9

6   0.000  0.000     80.3

7   0.000  0.000     80.0

8   0.000  0.000     79.7

9   0.000  0.000     79.8

10  1.414  0.000     78.4

11 -1.414  0.000     75.6

12  0.000  1.414     78.5

13  0.000 -1.414     77.2

• Applying Response Surface Method in R, using first-order model

> CR.rs1 <- rsm (response ~ FO(x1,x2), data=d1) 

> summary(CR.rs1)

Call:

rsm(formula = response ~ FO(x1, x2), data = d1)

            Estimate Std. Error  t value Pr(>|t|)    

(Intercept) 78.50000    0.37450 209.6107  < 2e-16 ***

x1           1.00755    0.47744   2.1103  0.06101 .  

x2           0.46735    0.47744   0.9789  0.35073    

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Multiple R-squared:  0.3511,  Adjusted R-squared:  0.2214 

F-statistic: 2.706 on 2 and 10 DF,  p-value: 0.115

Analysis of Variance Table

Response: response

            Df  Sum Sq Mean Sq F value    Pr(>F)

FO(x1, x2)   2  9.8671  4.9335  2.7058 0.1150132

Residuals   10 18.2329  1.8233                  

Lack of fit  6 18.0209  3.0035 56.6695 0.0007986

Pure error   4  0.2120  0.0530                  

Direction of steepest ascent (at radius 1):

       x1        x2 

0.9071625 0.4207806 

Corresponding increment in original units:

       x1        x2 

0.9071625 0.4207806 

Here the P-Value is 0.115. A p-value higher than 0.05 (> 0.05) is not statistically significant and indicates weak evidence against the null hypothesis. This means we fail to reject the null hypothesis and cannot accept the alternative hypothesis

• Applying Response Surface Method in R, using second-order model

> CR.rs2 <- rsm (response ~ SO(x1,x2), data=d1) 

> summary(CR.rs2)

Call:

rsm(formula = response ~ SO(x1, x2), data = d1)

             Estimate Std. Error  t value  Pr(>|t|)    

(Intercept) 79.939958   0.114754 696.6216 < 2.2e-16 ***

x1           1.007552   0.090728  11.1052 1.068e-05 ***

x2           0.467346   0.090728   5.1511  0.001322 ** 

x1:x2        0.225000   0.128299   1.7537  0.122924    

x1^2        -1.382706   0.097308 -14.2095 2.031e-06 ***

x2^2        -0.957578   0.097308  -9.8407 2.378e-05 ***

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Multiple R-squared:  0.9836,  Adjusted R-squared:  0.9719 

F-statistic: 83.96 on 5 and 7 DF,  p-value: 4.302e-06

Analysis of Variance Table

Response: response

            Df  Sum Sq Mean Sq  F value    Pr(>F)

FO(x1, x2)   2  9.8671  4.9335  74.9299 1.877e-05

TWI(x1, x2)  1  0.2025  0.2025   3.0755    0.1229

PQ(x1, x2)   2 17.5695  8.7848 133.4214 2.670e-06

Residuals    7  0.4609  0.0658                   

Lack of fit  3  0.2489  0.0830   1.5654    0.3293

Pure error   4  0.2120  0.0530                   

Stationary point of response surface:

       x1        x2 

0.3879028 0.2895971 

Eigenanalysis:

eigen() decomposition

$values

[1] -0.9296432 -1.4106413

$vectors

         [,1]       [,2]

x1 -0.2409913 -0.9705273

x2 -0.9705273  0.2409913