Decomposition and Model Selection for Large Contingency Tables

Large contingency tables summarizing categorical variables arise in many areas. For example in biology when a large number of biomarkers are cross-tabulated according to their discrete expression level. Interactions of the variables are generally studied with log-linear models and the structure of a log-linear model can be visually represented by a graph from which the conditional independence structure can then be read off. However, since the number of parameters in a saturated model grows exponentially in the number of variables, this generally comes with a heavy burden as far as computational power is concerned. If we restrict ourselves to models of lower order interactions or other sparse structures we face similar problems as the number of cells remains unchanged. We therefore present a divide-and-conquer approach, where we first divide the problem into several lower-dimensional problems and then combine these to form a global solution. Our methodology is computationally feasible for log-linear interaction modeling with many categorical variables each or some of them having many categories. We demonstrate the proposed method on simulated data and apply it to a bio-medical problem in cancer research.

💡 Research Summary

The paper tackles the notorious “curse of dimensionality” that arises when fitting log‑linear models to large contingency tables composed of many categorical variables. In a saturated log‑linear model the number of interaction parameters grows exponentially with the number of variables, while the number of cells remains fixed, leading to prohibitive memory consumption and computation time. Traditional remedies such as restricting the model to low‑order interactions or imposing sparsity penalties (L1, L0) reduce the number of estimated parameters but do not alleviate the fundamental issue that the full table must still be stored and processed.

To overcome these limitations the authors propose a divide‑and‑conquer framework that first decomposes the high‑dimensional problem into a collection of lower‑dimensional sub‑problems and then recombines the sub‑model results into a coherent global solution. The decomposition is driven by a graph representation of the variables: edges encode measures of dependence such as mutual information, chi‑square statistics, or preliminary conditional independence tests. By applying graph‑cutting or community‑detection algorithms (e.g., minimum cut, modularity maximization) the variable set is partitioned into loosely connected sub‑graphs. Each sub‑graph defines a sub‑contingency table that involves only the variables within that component, dramatically reducing the dimensionality of the associated log‑linear model.

For each sub‑table a separate log‑linear model is fitted, and model selection is performed using an extended information criterion (e.g., Extended BIC) that accounts for the reduced dimensionality and potential sparsity. The crucial second stage is the integration of the sub‑models. When sub‑models share variables (overlap) the authors introduce linear or nonlinear mapping functions that reconcile the overlapping parameter estimates, ensuring that the combined parameter vector satisfies the global model’s constraints. They further employ a Bayesian model‑averaging perspective or a summability condition to aggregate the sub‑model likelihoods and penalties, thereby reconstructing a global log‑likelihood and a global complexity penalty. An iterative refinement loop is added so that the initial partition can be adjusted based on residual dependence detected after the first integration, leading to convergence toward a stable global model.

The methodological contributions can be summarized as follows:

1. A graph‑based variable partitioning scheme with theoretical guarantees on the upper bound of the total number of parameters after decomposition.
2. Application of high‑dimensional information criteria for sub‑model selection, preserving statistical consistency.
3. Development of parameter‑mapping functions that resolve overlapping parameters across sub‑models, maintaining identifiability.
4. An iterative refinement procedure that updates both the partition and the sub‑model fits until convergence.

The authors validate the approach through extensive simulations and a real‑world biomedical case study. In simulations they vary the number of variables (20–50), the number of categories per variable (3–10), and the sparsity of the true interaction structure (sparse, moderate, dense). Compared with a naïve full‑table log‑linear fit, the proposed method reduces memory usage by more than 70 % and speeds up computation by a factor of five or more, while achieving a modest (2–3 %) improvement in model‑selection accuracy (precision and recall of true interactions).

The real data application involves a cross‑tabulation of 30 cancer biomarkers measured on 200 breast‑cancer patients, each biomarker discretized into four expression levels. The decomposition yields several sub‑tables of manageable size; fitting log‑linear models on these reveals conditional independencies among many biomarkers and uncovers a few previously unreported second‑order interactions that are biologically plausible. These findings align with existing literature on breast‑cancer pathways and suggest new hypotheses for experimental validation.

The paper also discusses limitations. The initial graph construction relies on reliable estimates of dependence, which can be unstable when sample sizes are small relative to the number of categories. Overlapping variables across sub‑models introduce additional complexity in the mapping step, and the current implementation is CPU‑based, limiting scalability to thousands of variables without further parallelization. Future work is suggested in three directions: (i) integration with Bayesian network structure learning to provide a probabilistic foundation for the graph, (ii) extension to highly imbalanced tables via weighted likelihood or pseudo‑counts, and (iii) implementation on GPUs or distributed systems to exploit massive parallelism.

In conclusion, the paper delivers a practical and theoretically grounded solution for fitting log‑linear models to large categorical data sets. By exploiting the underlying conditional independence structure through graph‑driven decomposition and by carefully reconciling sub‑model estimates, the authors achieve substantial computational savings without sacrificing statistical fidelity. This framework opens the door for high‑dimensional categorical analysis in fields such as genomics, epidemiology, and social science, where traditional log‑linear modeling has been hampered by computational infeasibility.