Simple Algorithm Picks Powerful Transcription Factor Teams for Gene Sets

Choosing which transcription factors to perturb in an experiment sounds straightforward until the arithmetic sets in. If a researcher wants a small team of factors that collectively governs a disease-related gene set, the possibilities multiply explosively: selecting just three factors from a pool of 400 already yields more than 10 million combinations. A new study in Algorithms for Molecular Biology tackles this combinatorial maze by recasting it as a stochastic probing problem on bipartite graphs, and reports that surprisingly simple, non-adaptive algorithms come within a few percent of the best achievable solution.

The work, by Lukas Geis, Dennis Hecker, Martin Hoefer, Ulrich Meyer and Marcel H. Schulz, addresses a gap in how regulatory networks are usually interrogated. Most existing approaches assess the importance of transcription factors one at a time, even though factors are known to act cooperatively — the classic example being the four Yamanaka factors that reprogram cell identity. Earlier efforts to design informative perturbation experiments, such as maximum-flow formulations or model-based tools like MEED, focus on maximizing knowledge about network relationships rather than on finding the factor set with the maximum coverage of regulated genes. The new framework, according to the authors, is the first to optimize the selection of a transcription factor subset explicitly for gene coverage.

The mathematical core is a weighted, complete bipartite graph in which one set of nodes represents regulators — the transcription factors — and the other represents positions, which in the biological application are genes. Each edge carries a weight describing how strongly a factor influences a gene, but crucially these weights are random variables with known distributions whose actual realizations are hidden until a regulator is “probed.” Probing a factor reveals the realized weights of all its incident edges at once. The task is to probe at most k regulators and then choose a subset of at most ℓ of them that maximizes an objective function evaluated on the realized weights.

Two broad classes of algorithms exist for such problems. A non-adaptive policy commits upfront to which k regulators it will probe, while an adaptive policy decides each next probe based on what previous probes revealed. Adaptive policies are formally stronger — the problem can be cast as a finite Markov decision process with an optimal, if computationally prohibitive, adaptive solution — but they are expensive to compute. The gap between the two is quantified by the adaptivity gap, the worst-case ratio of the expected objective value of an optimal adaptive policy to that of an optimal non-adaptive one. Building on a seminal result by Asadpour and Nazerzadeh, the authors note that when the relevant parent function is monotone submodular, this gap is bounded by roughly 1.58.

The paper explores three objective functions. Two of them, f_max and f_sum, treat regulators as fully independent, summing either the maximum or the total of a regulator’s incident edge weights; their parent functions are provably submodular, so the classical guarantees apply directly. The biologically realistic case is harder. The coverage function f_cov scores a chosen factor set by the total weight it covers, counting each gene only once under its strongest selected regulator — a generalization of the weighted Maximum Coverage problem. This function is itself submodular, but the authors prove that its parent function is not submodular whenever k exceeds ℓ plus one, demonstrated through a small scaled counterexample. That loss of structure matters: the coverage problem is NP-hard and cannot be approximated within a constant factor better than the standard threshold, even after all probes are complete.

Faced with these theoretical obstacles, the team designed modified greedy heuristics. For the adaptive side, AmpCov estimates the expected marginal gain of each unprobed regulator by treating unprobed edges at their expected values and running a greedy selection routine. The non-adaptive counterpart, NampCov, computes exact expected increases for the first ℓ rounds and then, when that becomes computationally intractable, falls back on picking the regulator with the highest expected sum of incident edge weights. General approximation factors for these modifications depend on the maximum node degree and, the authors concede, substantially underestimate what the algorithms actually achieve in practice.

The empirical evaluation ran on simulated networks under two regimes: a Uniform setting with equal numbers of regulators, positions and weight values, and a Network setting calibrated to resemble gene regulatory networks, with Poisson-distributed edge weights and many more genes than factors. For each configuration size, 100 random instances were generated and each algorithm was run ten times per parameter combination. In the Network setting with the coverage objective, NampCov achieved approximation factors between 0.94 and 0.99 relative to an offline optimum — itself an upper bound on what any adaptive policy could attain — nearly matching its adaptive rival everywhere except at k = ℓ, where a guarantee exists and the adaptive variant edged ahead. Larger instances and strategies that probed more than they kept both improved performance. Runtimes told a consistent story: for ℓ < k the non-adaptive algorithm was always faster, sometimes dramatically so, because the adaptive version must repeatedly invoke a greedy subroutine during each probing phase.

To demonstrate real-world utility, the researchers built transcription factor–gene networks from a healthy T-cell sample drawn from the International Human Epigenome Consortium EpiATLAS. Candidate regulatory elements came from H3K27ac ChIP-seq peaks, and element-to-gene links were predicted with an adapted Activity-by-Contact model using Hi-C contact data. Binding sites for 420 transcription factors — those with a ChIP-seq peak at their own promoter — were annotated with FIMO using motif collections from JASPAR, HOCOMOCO and other sources. An edge was assigned only when binding sites were enriched in a gene’s regulatory elements, judged by a Fisher’s exact test at a p-value threshold of 0.05, with the edge weight set to the fraction of the gene’s elements containing a binding site. Two gene sets were targeted: 131 genes annotated for T-cell mediated immunity and 83 genes linked to lymphoid leukemia from the DISGENET curated database.

Running the algorithms 1,500 times with the coverage objective, the team searched for transcription factor triplets with the best network coverage. The most frequently chosen triplets regulated far more genes of the networks than randomly assembled triplets, and the individual factors showed elevated specificity — they connected to a larger percentage of the target gene set than of the genome-wide background of annotated genes with at least one regulatory element. Co-expression among factors within the top triplets was significantly higher than in random triplets, tested on RNA-seq data from 1,555 EpiATLAS samples with a two-sided Mann–Whitney U test, hinting at shared regulatory programs. Literature checks bolstered the picks: ZBTB17, implicated in lymphoid development and T-cell phenotype specification; IRF1, required for the differentiation and function of specific T-cell subtypes; and PATZ1, important for multiple T-cell populations, dominated the immunity network. In the leukemia network, repeatedly chosen factors included ZNF354A, ZNF530, ZNF320, ZSCAN22 and the tumor suppressor KLF10, which also regulates T-cell subtype differentiation. Even so, coverage had limits: a set of ten factors captured only 53 percent of the immunity genes and 65 percent of the leukemia genes, reflecting genuine regulatory diversity within the gene sets.

The authors are candid about the caveats. Position weight matrices simplify binding by treating positions independently, binding sites were binarized, and the fraction-of-elements edge weight is unlikely to translate linearly into expression regulation, with cofactors, cooperative modes and chromatin modifications all unaccounted for. Genes are weighted equally, though importance-scaled weights could be added. The compensation is minimal data requirements: only a sample-specific measurement of regulatory element activity, such as ATAC-seq, DNase-seq or ChIP-seq of element-associated histone marks, is needed, and any bipartite network with positive edge weights can be substituted. The framework extends naturally to micro RNA–mRNA or ligand–receptor networks, and the code, along with the two T-cell networks, is available on GitHub. For perturbation study design, the method offers a computationally cheap, near-optimal starting point for choosing which transcription factors to test together.

Subject of Research: Stochastic probing algorithms for selecting transcription factor subsets that maximize gene coverage in regulatory networks

Article Title: Probing transcription factor subsets in gene regulatory networks

Article References: Geis, L., Hecker, D., Hoefer, M., Meyer, U., & Schulz, M. H. (2026). Probing transcription factor subsets in gene regulatory networks. Algorithms for Molecular Biology, 21(1), Article 7. https://doi.org/10.1186/s13015-026-00297-x

Image Credits: AI Generated

DOI: 10.1186/s13015-026-00297-x

Keywords: transcription factors, gene regulatory networks, bipartite graphs, stochastic probing, adaptive algorithms, adaptivity gap, submodular optimization, maximum coverage, T-cell immunity, lymphoid leukemia, H3K27ac ChIP-seq, computational biology