Reference Manuscript: Huang, Y., Li, K., Mandal, A. et al. ForLion: a new algorithm for D-optimal designs under general parametric statistical models with mixed factors. Stat Comput 34, 157 (2024).
This package implements the Forlion algorithm on the problem of designing an experiment with both discrete and continuous factors. The ForLion algorithm searches for locally optimal approximate designs under the D-criterion. When the true parameters are unknown, we can generate the EW D-optimal designs. This features allow users to create efficient and robust experimental designs that account for parameter uncertainty. The algorithm performs an exhaustive search in a design space with mixed factors while keeping high efficiency and reducing the number of distinct experimental settings. We mainly focus on generalized linear models (GLMs) and multinomial logit models (MLMs). Furthermore, implementing approximate designs with continuous settings can be challenging in practical applications. To address this issue, our package includes an algorithm that converts approximate designs into exact designs. This approximate-to-exact design algorithm ensures that the resulting exact design remains highly efficient while being more feasible for practical implementation.
This is an example of electrostatic discharge (ESD) experiment in Section 4 of the reference manuscript.
Lukemire et al. (2019) reconsidered the electrostatic discharge (ESD) experiment described by Whitman et al. (2006) with a binary response and five mixed factors. The first four factors LotA, LotB, ESD, Pulse take values in \(\{-1,1\}\), and the fifth factor Voltage \(\in [25, 45]\) is continuous.
The logistic regression model is:
\(\text{logit}(\mu)=\beta_0+\beta_1\text{LotA}+\beta_2\text{LotB}+\beta_3\text{ESD}+\beta_4\text{Pulse}+\beta_5\text{Voltage}+\beta_{34}(\text{ESD}\times\text{Pulse})\)
The parameter values: \(\boldsymbol \beta = (\beta_0, \beta_1, \beta_2, \beta_3, \beta_4, \beta_5, \beta_{34})=(-7.5, 1.50,-0.2,-0.15,0.25,0.35,0.4)\)
hfunc.temp = function(y) {c(y,y[4]*y[5],1);}; # y -> h(y)=(y1,y2,y3,y4,y5,y4*y5,1)
n.factor.temp = c(0, 2, 2, 2, 2) # 1 continuous factor with 4 discrete factors
factor.level.temp = list(c(25,45),c(-1,1),c(-1,1),c(-1,1),c(-1,1))
link.temp="logit"
beta.value = c(0.35,1.50,-0.2,-0.15,0.25,0.4,-7.5) # continuous first and intercept last to fit hfunc.tempThen, we use the ForLion_GLM_Optimal() function in ForLion package to search for the D-optimal design
ForLion_GLM_Optimal(n.factor=n.factor.temp, factor.level=factor.level.temp, hfunc=hfunc.temp,bvec=beta.value,link=link.temp,reltol=1e-8, rel.diff=0.03, maxit=500, random=FALSE, logscale=TRUE)
#> Design Output
#> ===============================================================
#> Count Design Point(s) Allocation
#> ---------------------------------------------------------------
#> 1 25.0000, -1.0000, -1.0000, -1.0000, -1.0000 0.0868
#> 2 25.0000, -1.0000, -1.0000, -1.0000, 1.0000 0.0460
#> 3 25.0000, -1.0000, -1.0000, 1.0000, -1.0000 0.1114
#> 4 25.0000, -1.0000, -1.0000, 1.0000, 1.0000 0.0856
#> 5 25.0000, -1.0000, 1.0000, -1.0000, -1.0000 0.0833
#> 6 25.0000, -1.0000, 1.0000, -1.0000, 1.0000 0.1005
#> 7 25.0000, -1.0000, 1.0000, 1.0000, -1.0000 0.0373
#> 8 25.0000, -1.0000, 1.0000, 1.0000, 1.0000 0.0923
#> 9 25.0000, 1.0000, -1.0000, 1.0000, -1.0000 0.0140
#> 10 25.0000, 1.0000, 1.0000, 1.0000, -1.0000 0.1328
#> 11 29.4156, -1.0000, -1.0000, -1.0000, 1.0000 0.0682
#> 12 29.0032, -1.0000, 1.0000, -1.0000, -1.0000 0.0128
#> 13 32.7255, -1.0000, 1.0000, 1.0000, -1.0000 0.0484
#> 14 32.6908, -1.0000, 1.0000, 1.0000, -1.0000 0.0763
#> 15 31.4224, -1.0000, -1.0000, 1.0000, -1.0000 0.0042
#> ===============================================================
#>
#> m:
#> [1] 15
#>
#> det:
#> [1] 1.264497e-05
#>
#> convergence:
#> [1] TRUE
#>
#> min.diff:
#> [1] 0.0347
#>
#> x.close:
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 32.7255 -1 1 1 -1
#> [2,] 32.6908 -1 1 1 -1
#>
#> itmax:
#> [1] 500This is an example of minimizing surface defects experiment, which is example in S.3 in the supplementary material of the reference manuscript.
A polysilicon deposition process for manufacturing very large-scale integrated (VLSI) circuits was described with six control factors, namely, Cleaning Method, Deposition temperature (\(^\circ\)C), Deposition pressure (mtorr), Nitrogen flow rate (seem), Silane flow rate (seem), and Settling time (minutes). Wu (2008) used the relevant experiment as an illustrative example and categorized the response Surface Defects from a count number to one of the five ordered categories, namely, Practically no surface defects (I, \(0\sim 3\)), Very few defects (II, \(4\sim 30\)), Some defects (III, \(31\sim 300\)), Many defects (IV, \(301\sim 1000\)), and Too many defects (V, \(1001\) and above).
The example uses cumulative proportional odds model.
\(\log\left(\frac{\pi_{i1}+\cdots+\pi_{ij}}{\pi_{i,j+1}+\cdots+\pi_{iJ}}\right)=\theta_{j} - \beta_1 x_{i1} - \beta_2 x_{i2} - \beta_3 x_{i3}\nonumber - \beta_4 x_{i4} - \beta_5 x_{i5} - \beta_6 x_{i6}\)
with \(i=1, \ldots, m\) and \(j=1, 2, 3, 4\) and \(\boldsymbol \theta = (\theta_1, \theta_2, \theta_3, \theta_4, \beta_1, \beta_2, \beta_3, \beta_4, \beta_5, \beta_6)=(-1.77994301, -0.05287782, 1.86852211, 2.76330779,\) $ -0.94437464, 0.18504420, -0.01638597, -0.03543202, -0.07060306, 0.10347917)$.
Note: Always put continuous factors ahead of discrete factors, pay attention to the order of coefficients pairing with predictors
link.temp = "cumulative"
## Note: Always put continuous factors ahead of discrete factors, pay attention to the order of coefficients paring with predictors
n.factor.temp = c(0,0,0,0,0,2) # 1 discrete factor w/ 2 levels + 5 continuous
factor.level.temp = list(c(-25,25), c(-200,200),c(-150,0),c(-100,0),c(0,16),c(-1,1))
J = 5 #num of response levels
p = 10 #num of parameters
hfunc.temp = function(y){
if(length(y) != 6){stop("Input should have length 6");}
model.mat = matrix(NA, nrow=5, ncol=10, byrow=TRUE)
model.mat[5,]=0
model.mat[1:4,1:4] = diag(4)
model.mat[1:4, 5] =((-1)*y[6])
model.mat[1:4, 6:10] = matrix(((-1)*y[1:5]), nrow=4, ncol=5, byrow=TRUE)
return(model.mat)
}
hprime.temp=NULL #use numerical gradient for optim, thus could be NULL, if use analytical gradient then define hprime function
b.temp = c(-1.77994301, -0.05287782, 1.86852211, 2.76330779, -0.94437464, 0.18504420, -0.01638597, -0.03543202, -0.07060306, 0.10347917)Then, we can use ForLion_MLM_Optimal() function in ForLion package to search for the D-optimal design
ForLion_MLM_Optimal(J=J, n.factor=n.factor.temp, factor.level=factor.level.temp, hfunc=hfunc.temp, h.prime=hprime.temp, bvec=b.temp, link=link.temp, Fi.func=Fi_MLM_func, delta=1e-2, epsilon=1e-12, reltol=1e-10, rel.diff=5e-4, maxit=500, optim_grad=FALSE)
#> Design Output
#> ==============================================================================
#> Count Design Point(s) Allocation
#> ------------------------------------------------------------------------------
#> 1 -25.0000, -200.0000, 0.0000, -18.1815, 16.0000, -1.0000 0.1237
#> 2 -25.0000, 200.0000, -150.0000, -59.9634, 16.0000, 1.0000 0.0689
#> 3 -25.0000, -200.0000, 0.0000, -14.8567, 0.0000, 1.0000 0.1132
#> 4 -22.7757, 200.0000, -150.0000, 0.0000, 0.0000, -1.0000 0.1208
#> 5 -25.0000, 200.0000, 0.0000, -86.0167, 16.0000, 1.0000 0.0915
#> 6 -25.0000, -51.3720, 0.0000, -100.0000, 0.0000, 1.0000 0.0269
#> 7 24.7015, 200.0000, 0.0000, 0.0000, 0.0000, -1.0000 0.0899
#> 8 -4.3339, 200.0000, 0.0000, -100.0000, 0.0000, 1.0000 0.0378
#> 9 -25.0000, 200.0000, -150.0000, -87.1986, 0.0000, 1.0000 0.0423
#> 10 -25.0000, 200.0000, 0.0000, -100.0000, 0.0000, -1.0000 0.0681
#> 11 -25.0000, 70.6646, 0.0000, -100.0000, 0.0000, -1.0000 0.0228
#> 12 -25.0000, 69.9565, 0.0000, -100.0000, 0.0000, -1.0000 0.0025
#> 13 -25.0000, -166.2424, -150.0000, 0.0000, 0.0000, 1.0000 0.0385
#> 14 -25.0000, 70.2402, 0.0000, -100.0000, 0.0000, -1.0000 0.0014
#> 15 -25.0000, -166.4595, -150.0000, 0.0000, 0.0000, 1.0000 0.0058
#> 16 -25.0000, 70.1238, 0.0000, -100.0000, 0.0000, -1.0000 0.0018
#> 17 -25.0000, -166.6505, -150.0000, 0.0000, 0.0000, 1.0000 0.0057
#> 18 -25.0000, 70.0653, 0.0000, -100.0000, 0.0000, -1.0000 0.0017
#> 19 -25.0000, -166.8095, -150.0000, 0.0000, 0.0000, 1.0000 0.0033
#> 20 -25.0000, -166.8887, -150.0000, 0.0000, 0.0000, 1.0000 0.0017
#> 21 -25.0000, -166.8881, -150.0000, 0.0000, 0.0000, 1.0000 0.0016
#> 22 -25.0000, -166.8871, -150.0000, 0.0000, -0.0000, 1.0000 0.0015
#> 23 -25.0000, -166.8851, -150.0000, 0.0000, 0.0000, 1.0000 0.0015
#> 24 -25.0000, -166.8832, -150.0000, 0.0000, 0.0000, 1.0000 0.0014
#> 25 -25.0000, -166.8825, -150.0000, 0.0000, 0.0000, 1.0000 0.0013
#> 26 -25.0000, -166.8790, -150.0000, 0.0000, 0.0000, 1.0000 0.0012
#> 27 -25.0000, -166.8785, -150.0000, 0.0000, 0.0000, 1.0000 0.0012
#> 28 -25.0000, -166.8744, -150.0000, 0.0000, 0.0000, 1.0000 0.0012
#> 29 -25.0000, -166.8703, -150.0000, 0.0000, 0.0000, 1.0000 0.0011
#> 30 -25.0000, -166.8639, -150.0000, 0.0000, 0.0000, 1.0000 0.0011
#> 31 -25.0000, -166.8630, -150.0000, 0.0000, 0.0000, 1.0000 0.0010
#> 32 -25.0000, 70.0295, 0.0000, -100.0000, 0.0000, -1.0000 0.0009
#> 33 -25.0000, -166.8529, -150.0000, 0.0000, 0.0000, 1.0000 0.0009
#> 34 -25.0000, -166.8493, -150.0000, 0.0000, 0.0000, 1.0000 0.0009
#> 35 -25.0000, -166.8434, -150.0000, 0.0000, 0.0000, 1.0000 0.0008
#> 36 -25.0000, 70.0376, 0.0000, -100.0000, 0.0000, -1.0000 0.0007
#> 37 -25.0000, -166.8382, -150.0000, 0.0000, 0.0000, 1.0000 0.0008
#> 38 -25.0000, -166.8270, -150.0000, 0.0000, 0.0000, 1.0000 0.0007
#> 39 -25.0000, -166.8218, -150.0000, 0.0000, 0.0000, 1.0000 0.0006
#> 40 -25.0000, -166.8083, -150.0000, 0.0000, 0.0000, 1.0000 0.0005
#> 41 -25.0000, 70.0543, 0.0000, -100.0000, 0.0000, -1.0000 0.0003
#> 42 -25.0000, -166.8046, -150.0000, 0.0000, 0.0000, 1.0000 0.0004
#> 43 -25.0000, -166.7956, -150.0000, 0.0000, 0.0000, 1.0000 0.0004
#> 44 -25.0000, 70.0608, 0.0000, -100.0000, -0.0000, -1.0000 0.0002
#> 45 -25.0000, 70.0615, 0.0000, -100.0000, 0.0000, -1.0000 0.0002
#> 46 -25.0000, -166.7908, -150.0000, 0.0000, 0.0000, 1.0000 0.0003
#> 47 -25.0000, 70.0624, 0.0000, -100.0000, -0.0000, -1.0000 0.0002
#> 48 -25.0000, -166.7871, -150.0000, 0.0000, 0.0000, 1.0000 0.0003
#> 49 -25.0000, 70.0637, 0.0000, -100.0000, -0.0000, -1.0000 0.0002
#> 50 -25.0000, 70.0642, 0.0000, -100.0000, -0.0000, -1.0000 0.0002
#> 51 -25.0000, 70.0648, 0.0000, -100.0000, 0.0000, -1.0000 0.0002
#> 52 -25.0000, -166.7840, -150.0000, 0.0000, 0.0000, 1.0000 0.0003
#> 53 -25.0000, -166.7816, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 54 -25.0000, 70.0661, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 55 -25.0000, -166.7794, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 56 -25.0000, 70.0678, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 57 -25.0000, -166.7789, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 58 -25.0000, -166.7783, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 59 -25.0000, 70.0686, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 60 -25.0000, -166.7765, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 61 -25.0000, -166.7759, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 62 -25.0000, 70.0694, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 63 -25.0000, -166.7750, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 64 -25.0000, 70.0705, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 65 -25.0000, 70.0717, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 66 -25.0000, -166.7705, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 67 -25.0000, -166.7691, -150.0000, 0.0000, 0.0000, 1.0000 0.0002
#> 68 -25.0000, 70.0730, 0.0000, -100.0000, -0.0000, -1.0000 0.0001
#> 69 -25.0000, -166.7664, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 70 -25.0000, -166.7648, -150.0000, 0.0000, -0.0000, 1.0000 0.0001
#> 71 -25.0000, 70.0753, 0.0000, -100.0000, -0.0000, -1.0000 0.0001
#> 72 -25.0000, -166.7629, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 73 -25.0000, 70.0763, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 74 -25.0000, -166.7618, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 75 -25.0000, -166.7606, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 76 -25.0000, 70.0776, 0.0000, -100.0000, 0.0000, -1.0000 0.0001
#> 77 -25.0000, -166.7594, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 78 -25.0000, -166.7581, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 79 -25.0000, -166.7574, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 80 -25.0000, -166.7561, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 81 -25.0000, -166.7556, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 82 -25.0000, 70.0810, 0.0000, -100.0000, -0.0000, -1.0000 0.0000
#> 83 -25.0000, -166.7543, -150.0000, 0.0000, -0.0000, 1.0000 0.0001
#> 84 -25.0000, -166.7534, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 85 -25.0000, -166.7527, -150.0000, 0.0000, 0.0000, 1.0000 0.0001
#> 86 -25.0000, 70.0825, 0.0000, -100.0000, 0.0000, -1.0000 0.0000
#> 87 -25.0000, -166.7512, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 88 -25.0000, -166.7505, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 89 -25.0000, -166.7499, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 90 -25.0000, 70.0830, 0.0000, -100.0000, -0.0000, -1.0000 0.0000
#> 91 -25.0000, 70.0838, 0.0000, -100.0000, 0.0000, -1.0000 0.0000
#> 92 -25.0000, -166.7478, -150.0000, 0.0000, -0.0000, 1.0000 0.0000
#> 93 -25.0000, 70.0847, 0.0000, -100.0000, 0.0000, -1.0000 0.0000
#> 94 -25.0000, -166.7460, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 95 -25.0000, 70.0853, 0.0000, -100.0000, 0.0000, -1.0000 0.0000
#> 96 -25.0000, -166.7448, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 97 -25.0000, -166.7441, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> 98 -25.0000, 70.0871, 0.0000, -100.0000, -0.0000, -1.0000 0.0000
#> 99 25.0000, 200.0000, 0.0000, 0.0000, 16.0000, 1.0000 0.1031
#> 100 -25.0000, -166.7424, -150.0000, 0.0000, 0.0000, 1.0000 0.0000
#> ==============================================================================
#>
#> m:
#> [1] 100
#>
#> det:
#> [1] 219352835
#>
#> convergence:
#> [1] FALSE
#>
#> min.diff:
#> [1] 5e-04
#>
#> x.close:
#> [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,] -25 -166.8790 -150 0 0 1
#> [2,] -25 -166.8785 -150 0 0 1
#>
#> itmax:
#> [1] 500