Generate central-composite designs and associated building blocks

These functions generate central-composite designs, or building blocks thereof. They allow for flexible choices of replications, aliasing of predictors and fractional blocks, and choices of axis or ‘star’ points.

Usage

cube(basis, generators, n0 = 4, reps = 1, coding, randomize = TRUE, 
    blockgen, bid = 1, inscribed = FALSE)
star(basis, n0 = 4, alpha = "orthogonal", reps = 1, randomize = TRUE)
dupe(design, randomize = TRUE, coding)
foldover(basis, variables, bid, randomize = TRUE)
ccd(basis, generators, blocks = "Block", n0 = 4, alpha = "orthogonal", 
    wbreps = 1, bbreps = 1, randomize = TRUE, inscribed = FALSE, 
    coding, oneblock = FALSE)

Arguments

basis

In cube and ccd, a formula, or an integer giving the number of variables. If the formula has a left-hand side, the variables named there are appended to the design and initialized to NA. In star, dupe, and foldover, basis is a coded.data object to use as a reference.

generators

Optional formula or list of formulas to generate aliased variables

n0

Integer giving the number of center points. In ccd, this can be a vector of two numbers for the numbers of center points in the cube blocks and the star blocks, respectively.

reps

Integer number of replications of the cube or the star. (This does not create replicate blocks; use djoin to do that.)

coding

List of coding formulas for the design variables (those in basis and generators). In dupe, coding may be used to change the coding formulas, e.g. in a situation where we want to use the same design as before but center it elsewhere.

randomize

Logical value determining whether or not to randomize the design. In ccd, each block is randomized separately.

blockgen

A formula, string, or list thereof. Each element is evaluated, and the distinct combinations define fractional blocks for the design. Unlike ccd, cube returns only one of these blocks.

bid

(For block ID.) An integer index (from 1 to number of blocks) of the fractional block to return. The indexes are defined by the standard ordering of the block generators; e.g. if blockgen is of length 2, the bid values of (1, 2, 3, 4) correspond to generated levels of (--, +-, -+, ++) respectively.

inscribed

Logical value; if FALSE, the cube points are at +/- 1 in each variable. If TRUE, the entire desgn is scaled down so that the axis points are at +/- 1 and the cube points are at interior positions. In cube only, inscribed may be given a numeric value: use the value of alpha anticipated when axis points are added; or use inscribed = TRUE to scale in anticipation of alpha = "spherical".

alpha

If numeric, the position of the ‘star’ points. May also be a character string that matches or partially matches one of these:

"orthogonal": the star points are positioned to block the design orthogonally
"rotatable": the star points are chosen to make the design rotatable
"spherical": the star points are the same distance as the corners of the design cube (alpha is the square root of the number of design factors)
"faces": the star points are face-centered (same as alpha = 1)

The user may specify a vector value of alpha if it is desired to vary them on different axes. The values are rotated cyclically as needed.

design

A coded.data object to be duplicated.

blocks

A string or a formula. If a character string, it is the name of the blocking factor; if a formula, the left-hand side is used as the name of the blocking factor. The formula(s) on the right-hand side are used to generate separate fractional blocks.

variables

Character vector of names of variables to fold over.

wbreps

Number(s) of within-block replications. If this is a vector of length 2, then separate numbers are used for the ‘cube’ and the ‘star’ blocks respectively.

bbreps

Number(s) of between-block replications (i.e., number of repeats of each block). If this is a vector of length 2, then separate numbers are used for the ‘cube’ and the ‘star’ blocks respectively.

oneblock

Logical. If TRUE, the blocking factor is removed and the whole design is randomized as a single block. Note that the default number of center points may be larger than you anticipated because they are combined.

Details

Central-composite designs (CCDs) are popular designs for use in response-surface exploration. They are blocked designs consisting of at least one ‘cube’ block (two-level factorial or fractional factorial, plus center points), and at least one ‘star’ block (points along each axis at positions -alpha and +alpha), plus center points. Everything is put on a coded scale, where the cube portion of the design has values of -1 and 1 for each variable, and the center points are 0.

The ccd function creates an entire CCD design; however, in practice, we often start with just the cube portion and build from there. Therefore, the functions cube, star, dupe, and foldover are provided, and one may use djoin to combine them.

In cube and ccd, the basis argument determines a basic design used to create cube blocks. For example, cube(basis = ~ A + B + C) would generate a basic design of 8 factorial points plus center points. Use generators if you want additional variables in a fractional design; for example, generators = c(D ~ -A*B, E ~ B*C) added to the above would generate a 5-factor design with defining relation I = -ABD = BCE = -ACDE. For convenience, basis may be an integer instead of a formula, in which case default variable names of x1, x2, ... are used; for example, cube(3, ~ -x1*x2*x3) generates a 1/2 fraction design with added center points.

If you want the cube points divided into fractional blocks, give the formula(s) in the blockgen argument of cube, or the blocks argument of ccd. For instance, suppose we call cube(basis = A+B+C+D+E, generators = F~-A*C*D). This design has 32 runs; but adding the argument blockgen = c("A*B*C","C*D*E") will create a fractional block of 32/4 = 8 runs. (cube is flexible; we could have used a formula instead, either blockgen = ~ c(A*B*C, C*D*E) or blockgen = c(~A*B*C, ~C*D*E).) Center points are added to each block as specified. In a call to ccd with the same basis and generators, adding blocks = Day ~ c(A*B*C, C*D*E) would do the same thing, only all 4 blocks will be included, and a factor named Day distinguishes the blocks.

The functions star, dupe, and foldover provide for creating new design blocks based on an existing design. They also provide for delayed evaluation: if the basis argument is missing, these functions simply return the call, djoin will fill-in basis = design1 and evaluate it.

dupe simply makes a copy of the design, and re-randomizes it. Therefore it is also a convenient way to re-randomize a design. If coding is provided, the coding formulas are replaced as well -- for example, to re-center the design.

Use star to generate star (axis) points, which consist of center points plus points at +/- alpha on each coordinate axis. You may specify the alpha you want, or a character argument to specify a certain criterion be met. For example, using delayed evaluation, ccd1 = djoin(cube1, star(alpha="sph")) will return a CCD with cube1 as the cube block, and with axis points at the same distance as the corners of the cube. Conditions for the criteria on alpha are described in detail in references such as Myers et al. (2009).

In star, determinations of orthogonality and rotatability are based on computed design moments of basis, rather than any assumptions about the structure of the design being augmented. Thus, it may be possible to augment an unusual design to obtain a rotatable design. Also, if an orthogonal star block is requested, the value of alpha may vary from axis to axis if that is required to satisfy the condition.

foldover reverses the levels of one or more design variables (i.e., those that are coded). By default, it reverses them all. However, if the bid argument is supplied, it instead returns the bidth fractional block that cube would have generated. That is, foldover(des, bid=3) is equivalent to cube(<arguments that created des>, bid=3) -- only it does so much more efficiently by folding on the appropriate factors.

In cases where there are constraints on the region of operability, you may want to specify inscribed = TRUE. This will scale-down the design so that no coded value exceeds 1. If using a building-block approach starting with a first-order design from cube, call cube with inscribed set to the anticipated value of alpha, or use inscribed = TRUE, and then use alpha = "spherical" in the subsequent call to star.

ccd generates an entire CCD. In practice, the building-block approach with cube, star, etc. is usually preferable, but ccd exists for convenience and backward compatibility with pre-2.00 versions of rsm. Many of the arguments are the same as those in cube; however, n0, wbreps, bbreps may be single values or vectors; if vectors, the first element is for the cube portions and the second element is for the star portions. In ccd, specifying wbreps is equivalent to specifying reps in a call to cube or star. bbreps refers to replicate blocks in the experiment, so that bbreps = c(2,3) specifies that we join two cube blocks and three blocks of star points.

If coding is not specified in a new design, default identity codings are created, e.g. x1 ~ x1.as.is.

Value

A coded.data object with the generated design, with additional variables run.order and std.order. If a multi-block design, the generated blocking variable will be a factor; all other variables will be numeric. The designs are sorted by blocks and run.order within blocks; and (unlike pre-1.41 versions of rsm) the row.names will be integers corresponding to this ordering. The user may sort by block and std.order within block to display the designs in their pre-randomized order.

References

Lenth RV (2009) ``Response-Surface Methods in R, Using rsm'', Journal of Statistical Software, 32(7), 1--17. doi:10.18637/jss.v032.i07

Myers, RH, Montgomery, DC, and Anderson-Cook, CM (2009) Response Surface Methodology (3rd ed.), Wiley.

Author

Russell V. Lenth

Note

Poor choices of generators and/or blocks can alias or partially alias some effects needed to estimate a second-order response surface. It is a good idea to run varfcn before collecting data to examine the prediction capabilities of the design and to ensure that the desired model can be fitted.

The function ccd.pick is available to help determine good choices for arguments to cube, star, and ccd.

An alternative to a CCD when you want to go straight to second-order modeling is a Box-Behnken design, generated by bbd. These designs are not as various or flexible as CCDs, but they can require fewer runs.

The non-exported function rsm:::.ccd.1.41 is provided in case it is needed by other packages for compatibility with old versions of rsm (version 1.41 or earlier). Given the same seed, it will also reproduce the randomization as a previously generated design from an old version.

Examples

library(rsm)

### Generate a standard 3-variable first-order design with 8 corner points and 4 center points
( FOdes <- cube (3, n0 = 4, coding = list (
                x1 ~ (Temp - 150)/10, x2 ~ (Pres - 50)/5, x3 ~ Feedrate - 4)) )
#>    run.order std.order Temp Pres Feedrate
#> 1          1        11  150   50        4
#> 2          2         5  140   45        5
#> 3          3         1  140   45        3
#> 4          4        10  150   50        4
#> 5          5        12  150   50        4
#> 6          6         9  150   50        4
#> 7          7         6  160   45        5
#> 8          8         7  140   55        5
#> 9          9         8  160   55        5
#> 10        10         3  140   55        3
#> 11        11         2  160   45        3
#> 12        12         4  160   55        3
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ (Temp - 150)/10
#> x2 ~ (Pres - 50)/5
#> x3 ~ Feedrate - 4
                                    
### Add an orthodonal star block with 12 runs to create a second-order CCD
( SOdes <- djoin(FOdes, star(n0=6)) )
#>    run.order std.order Temp Pres Feedrate Block
#> 1          1        11  150   50        4     1
#> 2          2         5  140   45        5     1
#> 3          3         1  140   45        3     1
#> 4          4        10  150   50        4     1
#> 5          5        12  150   50        4     1
#> 6          6         9  150   50        4     1
#> 7          7         6  160   45        5     1
#> 8          8         7  140   55        5     1
#> 9          9         8  160   55        5     1
#> 10        10         3  140   55        3     1
#> 11        11         2  160   45        3     1
#> 12        12         4  160   55        3     1
#> 13         1         6  150   50        6     2
#> 14         2         5  150   50        2     2
#> 15         3         3  150   40        4     2
#> 16         4         7  150   50        4     2
#> 17         5        11  150   50        4     2
#> 18         6         1  130   50        4     2
#> 19         7         8  150   50        4     2
#> 20         8        10  150   50        4     2
#> 21         9         4  150   60        4     2
#> 22        10         9  150   50        4     2
#> 23        11        12  150   50        4     2
#> 24        12         2  170   50        4     2
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ (Temp - 150)/10
#> x2 ~ (Pres - 50)/5
#> x3 ~ Feedrate - 4

### Same as above, except make the whole CCD at once; and make it rotatable
### and inscribed so that no coded value exceeds 1
SOdes2 <- ccd (3, n0 = c(4,6), alpha = "rotatable", inscribed = TRUE, coding = list (
                x1 ~ (Temp - 150)/10, x2 ~ (Pres - 50)/5, x3 ~ Feedrate - 4))

### Make two replicate blocks of FOdes (2nd one randomized differently)
djoin(FOdes, dupe(FOdes))
#>    run.order std.order Temp Pres Feedrate Block
#> 1          1        11  150   50        4     1
#> 2          2         5  140   45        5     1
#> 3          3         1  140   45        3     1
#> 4          4        10  150   50        4     1
#> 5          5        12  150   50        4     1
#> 6          6         9  150   50        4     1
#> 7          7         6  160   45        5     1
#> 8          8         7  140   55        5     1
#> 9          9         8  160   55        5     1
#> 10        10         3  140   55        3     1
#> 11        11         2  160   45        3     1
#> 12        12         4  160   55        3     1
#> 13         1         9  150   50        4     2
#> 14         2        11  150   50        4     2
#> 15         3         7  140   55        5     2
#> 16         4         1  140   45        3     2
#> 17         5         3  140   55        3     2
#> 18         6         2  160   45        3     2
#> 19         7         5  140   45        5     2
#> 20         8         6  160   45        5     2
#> 21         9         8  160   55        5     2
#> 22        10        10  150   50        4     2
#> 23        11         4  160   55        3     2
#> 24        12        12  150   50        4     2
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ (Temp - 150)/10
#> x2 ~ (Pres - 50)/5
#> x3 ~ Feedrate - 4

### Fractional blocking illustration (with no center points)
# Basic design (bid = 1 ---> block generators b1 = -1, b2 = -1)
block1 <- cube (~ x1 + x2 + x3 + x4,  generators = x5 ~ x1 * x2 * x3 * x4,
                n0 = 0, blockgen = ~ c(x1 * x2, x1 * x3), bid = 1)
block1                
#>   run.order std.order x1.as.is x2.as.is x3.as.is x4.as.is x5.as.is
#> 1         1         2       -1        1        1       -1        1
#> 2         2         1        1       -1       -1       -1       -1
#> 3         3         4       -1        1        1        1       -1
#> 4         4         3        1       -1       -1        1        1
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ x1.as.is
#> x2 ~ x2.as.is
#> x3 ~ x3.as.is
#> x4 ~ x4.as.is
#> x5 ~ x5.as.is

# The foldover (on all variables) of block1, in the same order
foldover(block1, randomize=FALSE)
#>   run.order std.order x1.as.is x2.as.is x3.as.is x4.as.is x5.as.is
#> 1         1         2        1       -1       -1        1       -1
#> 2         2         1       -1        1        1        1        1
#> 3         3         4        1       -1       -1       -1        1
#> 4         4         3       -1        1        1       -1       -1
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ x1.as.is
#> x2 ~ x2.as.is
#> x3 ~ x3.as.is
#> x4 ~ x4.as.is
#> x5 ~ x5.as.is

# The 4th fractional block:
( block4 <- foldover(block1, bid = 4) )
#>   run.order std.order x1.as.is x2.as.is x3.as.is x4.as.is x5.as.is
#> 1         1         2        1        1        1       -1       -1
#> 2         2         3       -1       -1       -1        1       -1
#> 3         3         4        1        1        1        1        1
#> 4         4         1       -1       -1       -1       -1        1
#> 
#> Data are stored in coded form using these coding formulas ...
#> x1 ~ x1.as.is
#> x2 ~ x2.as.is
#> x3 ~ x3.as.is
#> x4 ~ x4.as.is
#> x5 ~ x5.as.is