Choosing a design

Choosing a design#

Every design in this package answers a different question, and picking the wrong one is the most expensive mistake in the whole workflow — you find out only after the experiments are run. This page is the router: start from what you want to know, and it points at the design, the option settings that matter, and the notebook that works it through.

The design names below are exactly the ones in the Design type menu of the browser app, and the same names the CLI takes as discopt doe new --template NAME.

Start here#

What you want to know	Design	In the app
Which of these factors actually matter?	2-level factorial	`factorial-2level`
…and there are too many of them for every combination	Fractional factorial, Plackett-Burman	`fractional-factorial`, `plackett-burman`
…and some may act through curvature	Definitive screening	`definitive-screening`
…and some factors are hard to change, or runs come in batches	Split-plot, blocked factorial	Python: `split_plot_design`, `blocked_factorial_design`
…and something I cannot hold constant is also varying	Latin square family	`latin-square`, `graeco-latin`, `hyper-graeco-latin`
Where is the optimum, and what does the surface look like near it?	Central composite, Box-Behnken	`central-composite`, `box-behnken`
Nothing yet — I want coverage before I commit to a model	Latin hypercube	`latin-hypercube`
How precisely can I pin down this model’s coefficients?	Model-based optimal	`linear`, `polynomial-1d`, `response-surface-2d`, `response-surface-3d`
Same, but my factors are proportions of a blend	Scheffé mixture	`scheffe-linear`, `scheffe-quadratic`, `scheffe-special-cubic`
Same, but my model is not in that list	Symbolic	`symbolic`
What is the best operating point? (not: what is the model?)	Active learning	CLI: `discopt doe optimize`
Which of two rival models is right?	Discrimination design	Python: `discriminate_design`

Two rows are deliberately not in the browser app. Active learning needs a surrogate refit between rounds, and discrimination needs the base discopt solver; both are one command or one function call away at the keyboard. See the CLI workflow.

Which factors matter? — screening#

factorial-2level. Run every combination of each factor at a low and a high level. With \(k\) factors that is \(2^k\) runs, and from them you get a signed estimate of every main effect and every interaction, plus an ANOVA table [Box et al., 2005, Montgomery, 2017]. Levels may be numbers or strings, so “catalyst A vs catalyst B” is a factor like any other.

Use it first, before fitting any response surface. A quadratic model in five factors costs 21 coefficients; a screening run that shows three of the five do nothing turns that into a 10-coefficient model over a smaller box.

Centre points (numeric factors only) buy you two things for a handful of extra runs: a pure-error estimate that does not assume the model is right, and a curvature check — if the centre response sits well off the plane through the corners, a 2-level design cannot describe your system and you want a response surface.
Replicates repeat the whole design. Replicate when your measurement noise is the thing you are unsure about; add centre points when the model form is what you are unsure about.
\(k \gtrsim 7\) makes the full factorial unaffordable (128 runs). Drop to one of the designs below.

Worked through in Factor Screening with 2-Level Factorial Designs.

Many factors, few runs#

fractional-factorial runs a \(2^{k-p}\) fraction of the full factorial. Give the fraction by its generators, e.g. D=ABC, which works everywhere including the browser app. Or, in Python and at the command line, ask fractional_factorial_design to search for a fraction of a given size and resolution; the search needs the base discopt solver. Every fraction aliases some effects with others. alias_structure(design) prints the defining relation, the resolution and the alias groups, so check it before you run anything. fold_over(design) adds the mirror-image runs that break a resolution III design’s main-effect aliases.

plackett-burman estimates up to \(N-1\) main effects in \(N\) runs, where \(N\) is a multiple of 4 [Plackett and Burman, 1946]. Its main effects are only partially aliased with two-factor interactions, which spreads a real interaction thinly over many columns rather than onto one.

definitive-screening uses three levels and \(2m+1\) runs. Main effects are orthogonal to every two-factor interaction and every quadratic effect, so a curved factor is not mistaken for an inert one [Jones and Nachtsheim, 2011]. Analyse it with fit and the default linear model, then follow up on what is active.

A design with as many effects as runs leaves no residual degrees of freedom and so no ANOVA. discopt doe anova then reports the effects against Lenth’s pseudo standard error [Lenth, 1989], which assumes most effects are inert. lenth_pse and half_normal_scores give the same numbers in Python.

Something else is varying — blocking#

latin-square (3 factors), graeco-latin (4), hyper-graeco-latin (5). When a nuisance variable — day, operator, batch of feedstock, position in the furnace — will vary whether you like it or not, a square design arranges the runs so that its effect is estimated and removed rather than left in the residual [Cochran and Cox, 1957, Fisher, 1935]. One factor is the treatment you care about; the rest are blocks.

The cost is severe constraints on shape. Every factor must have the same number of levels \(k\); a Latin square is then \(k^2\) runs instead of \(k^3\), which is the whole point. Graeco-Latin squares need \(k \ge 3\) and \(k \ne 6\); hyper-Graeco-Latin needs three mutually orthogonal Latin squares, so \(k \in \{4, 5, 7, 8, \ldots\}\).

Warning

A single Latin square cannot test an interaction. Row and treatment together determine column, so the row×treatment subspace already contains the column main effect — the terms are aliased, and anova_report refuses the fit rather than reporting an F-ratio for a term that is partly some other term. Getting an interaction means giving something up: replicate the square and drop one block from the model, spending its degrees of freedom on the interaction. Replicating alone is not enough: with every block still in the model the interaction stays confounded (fully or partly, depending on how the second square was randomized) with a block, and anova_report either refuses it or warns that the sums of squares depend on term order. Shown in Latin-Square Designs and ANOVA.

Batches and hard-to-change factors. When a nuisance variable splits the runs of a factorial into batches, blocked_factorial_design confounds the highest-order interactions with the blocks, so the main effects stay clean. When a factor is slow or expensive to change, such as a furnace temperature, the runs are really a split-plot. The factor is set once per whole plot, and the other factors vary within it. split_plot_design builds that structure and split_plot_anova tests each effect against the right error term [Jones and Nachtsheim, 2009]. Analysing a split-plot as if it were fully randomized makes whole-plot effects look far more significant than they are.

Where is the optimum? — response surfaces#

Both designs here fit a full quadratic — intercept, main, square and cross terms — which is the cheapest model that has a stationary point to find [Box and Wilson, 1951, Montgomery, 2017].

central-composite is the classic: a \(2^k\) factorial core, \(2k\) axial points, and centre replicates.

Axial distance rotatable makes prediction variance depend only on distance from the centre, which is what you want when you do not know which direction matters. face puts the axial points on the face centres, giving a 3-level design.
By default every run stays inside the bounds you typed. Tick axial points outside bounds for the textbook scaling, where the axial runs sit beyond them — only if those settings are actually reachable.

box-behnken never visits a corner of the design box: every run sits at an edge midpoint, so no single run combines the extreme level of every factor [Box and Behnken, 1960]. Choose it when the corners are dangerous, infeasible, or simply outside what the equipment will do. It needs 3–5 factors (a central composite handles 2–6).

After the fit, canonical_analysis locates the stationary point and says whether it is a maximum, minimum, saddle or ridge. stationary_point_ci puts a confidence interval on its location [Del Castillo and Cahya, 2001]. A flat optimum performs well but is poorly located, and the interval says so. steepest_ascent_path gives the classical path out of a first-order region, and ridge_analysis gives the best point at each distance from the centre. To compare designs before running them, scaled_prediction_variance and fds_curve show how evenly each one predicts over the region [Zahran et al., 2003]. desirability and overall_desirability combine several responses into one [Derringer and Suich, 1980].

No model yet — space filling#

latin-hypercube. A stratified sample over the continuous box: each factor is split into \(n\) equal-probability strata and each stratum is used exactly once, so the sample covers every factor’s range no matter how many factors there are [McKay et al., 1979]. Unlike every other design here, you set the run count — it is not dictated by the factor count.

Reach for it when you want coverage without committing to a model form: fitting a surrogate, mapping a region before deciding where to look closely, or feeding an active-learning loop. The fit model option (linear or quadratic) only tells the later fit step which basis to estimate; it does not change the sample.

The space-filling criterion chooses among the many Latin hypercubes of the same size. discrepancy (the default) makes the points uniform. maximin pushes the closest pair of points as far apart as possible, which suits fitting a Gaussian-process surrogate [Morris and Mitchell, 1995]. space_filling_metrics scores any design on both criteria, and quasi_random_design gives Sobol and Halton sequences for comparison.

Estimate a known model precisely — optimal design#

linear, polynomial-1d, response-surface-2d, response-surface-3d.

These do not spread runs out for coverage — they place them where they carry the most information about the coefficients of a model you have already committed to, by maximizing a scalar function of the Fisher information matrix [Atkinson et al., 2007, Franceschini and Macchietto, 2008, Kiefer and Wolfowitz, 1959]. Expect the answer to look sparse and repetitive; a D-optimal design for a straight line puts half its runs at each end of the range, because that is genuinely the best way to estimate a slope.

Template	Model	Factors	Coefficients
`linear`	\(b_0 + \sum_i b_i x_i\)	1–8	\(1 + k\)
`polynomial-1d`	\(\sum_{j=0}^{d} b_j x^j\)	exactly 1	\(d + 1\)
`response-surface-2d`	full quadratic	exactly 2	6
`response-surface-3d`	full quadratic	exactly 3	10

Set the run count to at least the coefficient count — with fewer runs than coefficients the information matrix is singular and standard errors are meaningless. A little above it is where these designs pay off.

The theory is in Model-Based Design of Experiments.

Blends — mixture designs#

scheffe-linear, scheffe-quadratic, scheffe-special-cubic.

When the factors are proportions of a formulation they are not free: they sum to a fixed total, so the design space is a simplex rather than a box, and the usual polynomial has no intercept [Cornell, 2002, Scheffé, 1958, Scheffé, 1963]. Set mixture total to whatever your components sum to (1.0 for fractions, 100 for percentages).

With \(q\) components:

Template	Terms	Coefficients
`scheffe-linear`	pure blends only	\(q\)
`scheffe-quadratic`	+ pairwise blending	\(q + \binom{q}{2}\)
`scheffe-special-cubic`	+ three-way blending	\(q + \binom{q}{2} + \binom{q}{3}\)

Start linear unless you have reason to expect synergy between components; the quadratic terms are what a blend being better than either ingredient alone looks like. See Mixture Designs.

Component bounds. Real formulations bound each component, and the bounds carve an irregular region out of the simplex. check_mixture_bounds reports whether the bounds can be met at all and tightens the ones the others imply. extreme_vertices_design gives the classical design on the region: its vertices, edge and face centroids, and the overall centroid [McLean and Anderson, 1966]. to_pseudo_components rescales the region so the standard Scheffé designs apply again. The optimal-design templates accept the same bounds and refuse ones that no blend can satisfy.

Your own model#

symbolic. Type the response as an expression — k0 * exp(-Ea / (8.314 * T)) — list its parameters with nominal values, and it is differentiated symbolically to build a design that estimates those parameters as precisely as possible.

Important

A nonlinear model’s information depends on the parameter values, so the design is only optimal around the numbers you supply. This is not a flaw to work around, it is the loop: design around a guess, run it, fit, then re-centre the next batch on the fitted values with discopt doe extend [Franceschini and Macchietto, 2008, Galvanin et al., 2009]. Several turns of that loop are worked through in Model-Based Active Learning for Response Optimization.

Two practical notes. Parameter scale matters — an activation energy in J/mol is ~\(10^5\) while a pre-exponential may be ~1, and criteria computed from a matrix spanning ten orders of magnitude are numerically fragile; prefer parameters of comparable size, e.g. fitting \(E_a/R\) in kelvin. And check Identifiability and Estimability Analysis before designing for a model with many parameters: if two of them only ever appear as a product, no design will separate them, and the diagnostics say so before you spend the runs.

Beyond a single design#

Optimize, don’t model. If you want the best operating point and do not care about the model that gets you there, run active-learning rounds instead: fit a surrogate to what you have, propose the next batch by an acquisition function, repeat [Jones et al., 1998, Rasmussen and Williams, 2006, Snoek et al., 2012]. That is discopt doe optimize, worked through in Active-Learning Optimization with a Surrogate.

Several runs at once. With parallel capacity, design the batch jointly rather than one point at a time — the information matrix is additive over independent runs, so a joint design beats picking the best point \(N\) times [Galvanin et al., 2007, Sandrin et al., 2025]. See Batch / Parallel Design of Experiments.

Two rival models. When the question is which model rather than what parameters, design for the condition where the models disagree most relative to their prediction uncertainty [Buzzi-Ferraris et al., 1984, Hunter and Reiner, 1965]. See Model Discrimination Workflow.

Which criterion?#

The optimal-design templates take a criterion — a single number summarising “how much information”, since a matrix has no natural ordering [Atkinson et al., 2007].

Option	Optimizes	Choose it when
`determinant` (D)	maximize \(\log\det\mathrm{FIM}\)	Default. Minimizes the volume of the joint confidence ellipsoid — the best all-round answer, and invariant to rescaling the parameters.
`trace` (A)	minimize \(\operatorname{tr}\mathrm{FIM}^{-1}\)	You care about the average of the individual parameter variances rather than their joint volume.
`min_eigenvalue` (E)	maximize \(\lambda_{\min}\)	One combination of parameters is much worse determined than the rest and you want to fix the worst case.
`condition_number` (ME)	minimize \(\lambda_{\max}/\lambda_{\min}\)	The parameters are badly correlated and you want a balanced design more than a maximally informative one.

Start with determinant. Switch to min_eigenvalue or condition_number when a fit comes back with one enormous standard error, or with two parameters whose errors are nearly perfectly correlated — that is the geometry those criteria are for.

How many runs?#

Screening: \(2^k\), times replicates, plus centre points. Non-negotiable — the structure is the design.
Latin square family: \(k^2\) per replicate, where \(k\) is the shared level count.
Response surface: set by the template (\(2^k + 2k +\) centres for a central composite).
Latin hypercube: your call. A common starting point for surrogate work is 10 runs per factor.
Optimal designs: at least the coefficient count from the tables above, and more if you want to test lack of fit rather than only estimate coefficients.

Everywhere: below the parameter count you get no standard errors, and at exactly the parameter count you get a perfect fit with zero residual degrees of freedom, which tells you nothing about whether the model is right.

After the experiments#

Fill in the response column and drop the workbook back on the app — the fit and the ANOVA run on their own. Both are reported where they apply: a model-based design gets a coefficient table with confidence intervals, a Latin square gets the factor-level F-table, and whichever does not apply says why. At the command line the same steps are discopt doe fit and discopt doe anova, plus discopt doe extend to design the next batch from what you now know (Design of Experiments from the Command Line).

Four ways this goes wrong#

Reading an interaction out of a single Latin square. It is aliased with a block; the guard refuses it. Replicate and give up a block.
Trusting a locally-optimal design far from the truth. For a nonlinear model the design is optimal around your nominal parameters and nowhere else. Design, fit, re-centre — do not run a hundred points off one guess.
Fitting different conditions as replicates. In a sequential loop, a runner that always returns the response under the same key hands the estimator a pile of values at one nominal condition, throwing away the very conditions the design chose. Give each observation its own response key (Batch / Parallel Design of Experiments).
Skipping the identifiability check. Diagnostics on the model you already have cost nothing; a batch of experiments that cannot separate two parameters costs a week (Identifiability and Estimability Analysis).