```{r} library(rsm) library(daewr) ``` ## Setting up a CCD with `ccd` ![CCD Summary](figures/ccd_summary.png) ```{r} ccd(y ~ x1 + x2, n0=c(5,0), alpha="rotatable", randomize=FALSE, oneblock=TRUE) ``` A CCD can be split into blocks. The most common split is two blocks: * Block 1: Factorial points + center points * Block 2: Axial points + center points ```{r} ccd(y ~ x1 + x2, n0=c(4,4), alpha="rotatable", randomize=FALSE, oneblock=FALSE) ``` For uniform precision, the total number of center points must equal or exceed the recommend number. In practice, the center points are split equally, with a few extra center points added to both blocks. It is also possible to split a CCD into three blocks. All the axial points remain in one block, and the factorial core is split into two. Center points are added to each block. ## Example: Optimizing the yield of a chemical reaction Time varies from 80 - 90 minutes (-1 to 1) Temp varies from 170 - 180 C (-1 to 1) The coded variable $$x$$ can be calculated using the original variable $$\xi$$ and the original values at +1 ($$\xi_{+1}$$) and -1 ($$\xi_{-1}$$). $$ x = \frac{\xi - \mathrm{mean}(\xi_{+1},\xi_{-1})}{(\xi_{+1} - \xi_{-1})/2} $$ ```{r} react <- ccd(Yield ~ x1 + x2, n0 = c(3,3), randomize=FALSE, alpha="rotatable", coding = list(x1 ~ (Time - 85)/5, x2 ~ (Temp - 175)/5)) react ``` ```{r} varfcn(react, ~ SO(x1,x2)) ``` ```{r} varfcn(react, ~ SO(x1,x2), contour=TRUE) ``` Let's add the measured yield data. ```{r} react$Yield <- c(80.5, 82.0, 81.5, 83.5, 83.9, 84.3, 84.0, 75.6, 78.4, 77.0, 78.5, 79.7, 79.8, 79.5) react ``` Fitting a FO+TWI model using the first 7 runs (FF + 3 CPs) ```{r} react_first <- react[1:7, ] FOreact <- rsm(Yield ~ FO(x1,x2) + TWI(x1,x2), data=react_first) summary(FOreact) ``` ```{r} # SO(x1,x2) = FO(x1,x2) + TWI(x1,x2) + PQ(x1,x2) SOreact <- rsm(Yield ~ Block + SO(x1,x2), data=react) summary(SOreact) xs(SOreact) ``` ## Visualizing the results ```{r} contour(SOreact, x1 ~ x2, image=TRUE) ``` ```{r} persp(SOreact, x1 ~ x2, contour=TRUE, zlab="Yield") ```