I recently became the manager of a complex piece of laboratory equipment. It has many knobs and dials to turn, and the different knobs and dials affect the significantly affect the responses I try to measure. Since there was no theory of how exactly each knob affects the response, I turned to statistical design of experiments figure out how to set my knobs. Using the same exact input each time, I use different settings for the various knobs and am trying to find the optimal (in my case, the maximal) setting for the knobs given this identical, unvarying input.
I initially expected that many settings would have an optimum level somewhere within the design space I was analyzing. Thus, I decided to use a definitive screening design to design my runs. I did exactly this, measuring the same sample nineteen different ways, with each of the nineteen runs having one of the nine knobs set to either a "high", "medium", or "low" setting. This was straightforward, using the the DefScreen()
function from the daewr
package in R
. My (normalized) responses as well as the (scaled and centered) experimental conditions are displayed below.
>>> scaled_data %>% dput
structure(list(method_num = 1:19, y = c(0.90085050161025, 0.00973935039938332,
0.046907811292605, 0.0477583584054755, 0.000734906456458273,
1, 0.101817637536817, 0.140301184447718, 0, 0.421559489775223,
0.636661367259339, 0.00451827441824333, 0.0714993448043543, 0.0605393859587476,
0.115095370854818, 0.222299049967782, 0.176960666274614, 0.0455311915787348,
0.203887872677827), A = c(0, 0, 1, -1, -1, 1, -1, 1, 1, -1, -1,
1, 1, -1, -1, 1, -1, 1, 0), B = c(1, -1, 0, 0, 1, -1, -1, 1,
-1, 1, -1, 1, 1, -1, -1, 1, 1, -1, 0), C = c(1, -1, 1, -1, 0,
0, 1, -1, 1, -1, -1, 1, -1, 1, -1, 1, 1, -1, 0), D = c(1, -1,
-1, 1, -1, 1, 0, 0, -1, 1, -1, 1, -1, 1, 1, -1, -1, 1, 0), E = c(1,
-1, 1, -1, 1, -1, 1, -1, 0, 0, 1, -1, 1, -1, 1, -1, -1, 1, 0),
F = c(1, -1, -1, 1, -1, 1, -1, 1, 1, -1, 0, 0, 1, -1, 1,
-1, 1, -1, 0), G = c(1, -1, -1, 1, 1, -1, -1, 1, 1, -1, 1,
-1, 0, 0, -1, 1, -1, 1, 0), H = c(1, -1, 1, -1, -1, 1, -1,
1, -1, 1, 1, -1, -1, 1, 0, 0, 1, -1, 0), J = c(1, -1, -1,
1, -1, 1, 1, -1, -1, 1, 1, -1, 1, -1, -1, 1, 0, 0, 0)), class = "data.frame", row.names = c(NA,
-19L), .Names = c("method_num", "y", "A", "B", "C", "D", "E",
"F", "G", "H", "J"))
Some ad-hoc model-fitting of this data revealed a few surprises:
- First and foremost, no quadratic terms seemed significant.
- Three (possibly four) of the variables of my original nine seemed the most important.
- One model that fits my data well is
y ~ H*J + B + B:G + B:G:H
- Even with those interactions included, there isn't clear evidence for any distinct optima within my design space.
A graph of my data vs. model predictions is below.
My main question is, what should I do next? Practically speaking, it isn't much work or time for me to do further runs, and so agonizing over the "best" way to proceed from the point of view of statistics not necessarily worth tons of time...yet I am curious about statistics and design of experiments in general. I figured the CV community would be a good place to ask for opinions and guidance on "official" ways to proceed. So, suppose for the sake of argument that further experimentation was very expensive and time-consuming. Some specific questions I have include:
The modeling I have done so far suggests (to me) that my design space didn't span any strong (local) optima. Thus, I need to extrapolate. How should I choose to extrapolate my original design space?
Should I use D-optimal design methods to augment my initial design?
How (if at all) should I use my model from my runs so far to inform future runs?
Can I just decide that the four variables I've identified as significant are all that matters, and re-do a definitive screening design centered at the point I identified as optimal in my first round of experiments?
What am I missing that I should be thinking about? My goal is to design an efficient set of experiments that will get me to the "optimum" setting for the knobs and dials on my experimental machine.