I need to optimize a production process by way of adjusting a total of 6 factors, all of them numerical variables with 3 levels. I expect the responses to be related to both the main effects and their quadratic terms, and I expect some of the second order interactions (3 of them, to be precise) to have an effect as well.
Based on this, I figured out I would use AlgDesign and the OptFederov function to get a D-optimal design.
I created the corresponding full factorial design matrix and formulated the model I want to attempt to fit to the data:
df <- gen.factorial(levels = 3, nVars = 6, varNames = LETTERS[1:6])
mod <- formula('~ A + A^2 + B + B^2 + C + C^2 + D + D^2 + E + E^2 + F + F^2 + + A:C + A:D + C:E + C:F + E*F')
Using OptFedorov I get:
design1 <- optFederov(frml = mod, data = df, nullify = 1)
unique(design1$design$A)
[1] -1 1
The same thing applies to all variables, that is, OptFedorov seems to only consider the extreme levels, not the mid point of each factor which, as far as I understand it, would not allow me to fit a quadratic function at all.
Since the documentation for AlgDesign suggests to use the quad() function for models with quadratic terms, I tried just that, and:
mod2 <- formula('~A + B + C + D + E + F + + A:C + A:D + C:E + C:F + E:F')
design2 <- optFederov(frml = quad(mod2), data = df, nullify = 1)
Error in quad(mod2) : could not find function "quad"
Using expand.formula() doesn't help either:
expand.formula(~.2, colnames(variables2))
~0.2
All the documentation I read makes use of the quad() function or some variation of the dot notation in the formula. So, am I missing something with the formulation of the model, or am I getting lost in localization or incompatible versions of R?
For reference the output of sessionInfo() is>
R version 4.2.0 (2022-04-22 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)
Matrix products: default
Random number generation:
RNG: Super-Duper
Normal: Inversion
Sample: Rejection
locale:
[1] LC_COLLATE=Italian_Italy.utf8 LC_CTYPE=Italian_Italy.utf8 LC_MONETARY=Italian_Italy.utf8
[4] LC_NUMERIC=C LC_TIME=Italian_Italy.utf8
attached base packages:
[1] grid stats graphics grDevices utils datasets methods base
other attached packages:
[1] tibble_3.1.7 AlgDesign_1.2.1 agricolae_1.3-5 FrF2_2.3 DoE.base_1.2-2
[6] conf.design_2.0.0
loaded via a namespace (and not attached):
[1] Rcpp_1.0.8.3 vcd_1.4-11 lattice_0.20-45 zoo_1.8-10 digest_0.6.29
[6] lmtest_0.9-40 utf8_1.2.2 gmp_0.7-1 mime_0.12 R6_2.5.1
[11] labelled_2.10.0 evaluate_0.15 highr_0.9 pillar_1.7.0 Rdpack_2.4
[16] rlang_1.0.2 rstudioapi_0.13 minqa_1.2.4 miniUI_0.1.1.1 nloptr_2.0.3
[21] Matrix_1.4-1 combinat_0.0-8 mathjaxr_1.6-0 splines_4.2.0 partitions_1.10-7
[26] lme4_1.1-29 polynom_1.4-1 stringr_1.4.0 questionr_0.7.7 numbers_0.8-5
[31] igraph_1.4.2 shiny_1.7.1 compiler_4.2.0 httpuv_1.6.5 xfun_0.31
[36] pkgconfig_2.0.3 htmltools_0.5.2 tidyselect_1.1.2 fansi_1.0.3 crayon_1.5.1
[41] dplyr_1.0.9 later_1.3.0 MASS_7.3-56 rbibutils_2.2.13 nlme_3.1-157
[46] xtable_1.8-4 lifecycle_1.0.1 magrittr_2.0.1 stringi_1.7.3 cli_3.3.0
[51] promises_1.2.0.1 scatterplot3d_0.3-41 ellipsis_0.3.2 generics_0.1.2 vctrs_0.4.1
[56] boot_1.3-28 klaR_1.7-1 tools_4.2.0 forcats_0.5.1 glue_1.6.2
[61] purrr_0.3.4 sfsmisc_1.1-15 hms_1.1.1 fastmap_1.1.0 colorspace_2.0-3
[66] cluster_2.1.3 knitr_1.39 haven_2.5.0