> ## Ranitidine Experiment, table 9.12
> dat=read.table("table9_12.dat", h=F)
> dat
   V1    V2    V3    V4
1   1  0.00 -1.41 6.943
2   2  1.00 -1.00 6.248
3   3 -1.41  0.00 2.100
4   4  0.00  0.00 2.034
5   5  0.00  0.00 2.009
6   6  0.00  0.00 2.022
7   7  1.00  1.00 3.252
8   8  1.41  0.00 9.445
9   9  0.00  1.41 1.781
10 10  0.00  0.00 1.925
11 11 -1.00 -1.00 2.390
12 12 -1.00  1.00 2.066
13 13  0.00  0.00 2.113
> A=dat[,2]
> B=dat[,3]
> lnCEF=dat[,4]
> g=lm(lnCEF ~ (A+B)^2 +I(A^2) +I(B^2))  # model (9.5) using all runs 
> summary(g)  # Table 9.13

Call:
lm(formula = lnCEF ~ (A + B)^2 + I(A^2) + I(B^2))

Residuals:
      Min        1Q    Median        3Q       Max 
-1.962024 -0.151854 -0.002400  0.090763  1.748289 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)   
(Intercept)   2.0244     0.5524   3.665  0.00802 **
A             1.9308     0.4374   4.414  0.00310 **
B            -1.3288     0.4374  -3.038  0.01890 * 
I(A^2)        1.4838     0.4703   3.155  0.01605 * 
I(B^2)        0.7743     0.4703   1.646  0.14371   
A:B          -0.6680     0.6176  -1.082  0.31529   
---
Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 

Residual standard error: 1.235 on 7 degrees of freedom
Multiple R-Squared: 0.8554,	Adjusted R-squared: 0.7521 
F-statistic: 8.282 on 5 and 7 DF,  p-value: 0.007461 

> ## Fig. 9.10 estimated response surface
> x=seq(-2, 2, .01)
> y=seq(-2, 2, .01)
> f=function(x, y)  2.0244+ 1.9308*x -1.3288*y + 1.4838*x^2 + 0.7743*y^2 -0.6680*x*y   ## (9.21)
> z=outer(x, y, f)
> contour(x*(5.25-4.5)+5.25, y*(20-14)+20, z)   # Fig. 9.10