Markov basis for design of experiments with three-level factors

We consider Markov basis arising from fractional factorial designs with three-level factors. Once we have a Markov basis, $p$ values for various conditional tests are estimated by the Markov chain Monte Carlo procedure. For designed experiments with a single count observation for each run, we formulate a generalized linear model and consider a sample space with the same sufficient statistics to the observed data. Each model is characterized by a covariate matrix, which is constructed from the main and the interaction effects we intend to measure. We investigate fractional factorial designs with $3^{p-q}$ runs noting correspondences to the models for $3^{p-q}$ contingency tables.

💡 Research Summary

The paper addresses a fundamental problem in the analysis of designed experiments where the response in each experimental run is a single count, such as the number of defects, events, or occurrences. Traditional analysis of variance techniques assume continuous responses and rely on asymptotic chi‑square approximations that become unreliable when the data are sparse or when the model includes high‑order interactions. To overcome these limitations, the authors develop a rigorous framework that combines algebraic statistics—specifically Markov bases—with fractional factorial designs involving three‑level factors.

The authors begin by formulating the experimental data within a generalized linear model (GLM) framework. The response vector (y) follows a Poisson (or binomial) distribution, and the design matrix (X) encodes the main effects and selected interaction effects of the three‑level factors. The sufficient statistics for the model are the marginal totals (t = X^{\top} y). The conditional distribution of (y) given (t) is the target for exact hypothesis testing, but its normalising constant is intractable for all but the simplest designs.

To sample from this conditional distribution, the paper leverages the concept of a Markov basis. A Markov basis is a finite set of integer moves ({d^{(1)},\dots,d^{(K)}}) such that, for any two tables (y) and (y’) sharing the same sufficient statistics, there exists a sequence of moves that transforms (y) into (y’) while never leaving the non‑negative orthant. The existence of such a basis guarantees that a simple Metropolis‑Hastings chain, which proposes adding or subtracting a randomly chosen move, is irreducible and aperiodic on the fiber ({y: X^{\top}y = t}). Consequently, the chain converges to the exact conditional distribution, allowing unbiased estimation of p‑values for any test statistic.

A key contribution of the paper is the identification of a structural equivalence between fractional factorial designs with (3^{p-q}) runs and (3^{p-q}) contingency tables. By coding the three levels of each factor as ({0,1,2}) and interpreting each experimental run as a cell in a multi‑dimensional table, the authors show that the marginal constraints imposed by the design matrix correspond precisely to the usual marginal constraints of a contingency table. This insight permits the direct application of existing results from algebraic statistics on toric ideals and Gröbner bases to the experimental design context. In particular, the toric ideal generated by the columns of (X) captures all algebraic relations among the cell counts, and its minimal generating set yields a Markov basis.

The computational pipeline described in the paper proceeds as follows: (1) specify the model (i.e., which main effects and interactions are of interest) and construct the corresponding design matrix (X); (2) compute the observed sufficient statistics (t); (3) use software such as 4ti2 to compute a Gröbner basis for the toric ideal (I_X), which provides the Markov basis; (4) run a Metropolis‑Hastings chain that proposes moves from the basis, accepting only those that keep all cell counts non‑negative; (5) after a burn‑in period, collect a large number of samples to approximate the conditional distribution of the chosen test statistic; (6) estimate the exact p‑value as the proportion of sampled statistics at least as extreme as the observed one.

Two illustrative examples are presented. The first involves a (3^{4-2}=81)-run design with four factors, where the model includes all main effects and all two‑factor interactions. The resulting Markov basis contains 27 moves, and the MCMC sampler mixes rapidly, yielding a precise p‑value that differs noticeably from the chi‑square approximation due to several zero‑count cells. The second example examines a larger (3^{5-3}=243)-run design with five factors and a model that includes selected three‑factor interactions. Here the Markov basis expands to 112 moves, yet the authors demonstrate that the chain still converges within a reasonable computational budget, confirming the feasibility of the approach for moderately large designs.

The discussion highlights several important practical considerations. First, the size of the Markov basis grows quickly with the number and order of interactions, which can increase memory usage and slow down the MCMC. The authors therefore recommend a careful selection of interaction terms based on scientific relevance, thereby keeping the basis manageable. Second, while the paper focuses on three‑level factors, the algebraic framework extends naturally to higher‑level factors, though the associated toric ideals become more complex. Third, for very large designs, approximate methods such as using a subset of “primitive” moves or employing sequential importance sampling may be necessary.

In conclusion, the paper provides a comprehensive methodological bridge between algebraic statistics and the design of experiments with discrete responses. By establishing the equivalence between three‑level fractional factorial designs and multi‑dimensional contingency tables, it enables the construction of exact Markov bases, which in turn facilitate rigorous conditional testing via MCMC. This approach offers a powerful alternative to asymptotic methods, especially in settings with sparse data or complex interaction structures, and opens new avenues for exact inference in modern experimental science.