Symbolic Regression Overview

• Symbolic regression is a method that discovers explicit mathematical expressions by jointly optimizing data fit and model complexity.
• It leverages evolutionary strategies like genetic programming and deep learning techniques such as transformer-based models to navigate vast search spaces.
• This approach uncovers underlying scientific laws in fields like physics, biology, and engineering, offering interpretable and actionable insights.

Symbolic regression (SR) is the task of discovering explicit mathematical expressions that relate input variables to measured outputs, using only data rather than specifying the form a priori. Unlike classical regression, which prescribes a fixed function class (e.g., polynomials of bounded degree, exponentials, or neural networks), SR jointly searches over both the structure of symbolic expressions and their numerical parameters. The overarching goal is to identify compact, interpretable models that generalize and often reveal underlying scientific laws or phenomenological relationships.

1. Formalism and Problem Setting

Symbolic regression seeks, for a dataset $\mathcal{D} = \{(\mathbf{x}_i, y_i)\}_{i=1}^N$ with $\mathbf{x}_i\in\mathbb{R}^d$, $y_i\in\mathbb{R}$, a symbolic expression $f(\mathbf{x})$ in a function space $\mathcal{S}$ (compositions of operators on $\mathbf{x}$ and constants) that minimizes a loss function measuring data fit, commonly mean squared error:

$$f^* = \arg\min_{f\in\mathcal{S}} \sum_{i=1}^N (f(\mathbf{x}_i) - y_i)^2 + \lambda\,\mathrm{Complexity}(f)$$

The complexity term—frequently the size or depth of the expression tree, or a description length in bits—encourages interpretable, parsimonious solutions (Makke et al., 2022, Bartlett et al., 17 Dec 2025). Expressions are typically represented as rooted, labeled trees: internal nodes are operators ($+$, $-$, $\sin$, $\exp$, etc.); leaves are variables $x_j$ or real constants.

This optimization constitutes a combinatorial mixed discrete–continuous problem known to be NP-hard—even under severe restrictions on the operator set and loss metric (Virgolin et al., 2022). For example, even for $\mathcal{S}$ comprising solely sums of variables (i.e., $f(\mathbf{x}) = \sum_k c_k x_k$ with integer $c_k \geq 0$), the decision version (can error be driven to zero?) is NP-complete via reduction from Unbounded Subset Sum.

2. Classical Algorithms: Genetic Programming and Its Variants

The canonical solution paradigm for SR is genetic programming (GP), which maintains a population of random expression trees. These are iteratively evolved using variation operators:

• Crossover: Exchange of subtrees between parent trees.
• Mutation: Stochastic replacement or modification of subtrees.
• Selection: Individuals are chosen for reproduction based on fitness (data fit penalized by complexity).
• Elitism: High-fitness solutions are retained across generations.

Symbolic expressions are decoded by tree traversal. GP offers generality—combinatorially searching over arbitrary compositions—but it is non-deterministic, can require many function evaluations (high computational cost), and is susceptible to “bloat” (the proliferation of large, functionally redundant trees) (Makke et al., 2022, Wang et al., 2019).

Modern variants of GP address some of these deficits:

• Geometric Semantic GP: Operators act in output (“semantic”) space, promoting smoother search landscapes (Wang et al., 2019).
• Multiple Regression GP (MRGP): Pools subtrees for subsequent selection via regression.
• Cartesian GP: Encodes expressions as acyclic graphs, improving subexpression reuse.
• Linear/Polynomial/Sparse GP: Constrains the model class for tractability or solution interpretability.

GP is augmented or hybridized with numerical solvers for efficient real constant fitting (e.g., Levenberg–Marquardt), Bayesian kernels (RVM), or deterministic sparse regression (FFX).

Uniform random global search (SRURGS) offers a baseline: it samples expressions uniformly from the grammar up to a size bound and optimizes constants, serving as a control to judge whether more sophisticated heuristics add value (Towfighi, 2019). While SRGP nearly always outperforms SRURGS on simple grammars, SRURGS can be more robust on complex, heterogeneous operator sets due to its lack of inductive bias.

3. Deep Learning–Based and Generative Methods

Recent developments leverage deep generative models and sequence-to-sequence architectures to scale SR:

• Transformer-based models (e.g., SymbolicGPT): Treat expression recovery as conditional language modeling: $p(e \mid D) = \prod_t p(e_t \mid e_{<t}, D)$. The data $D$ are embedded (e.g., via a T-net for order-invariance), and a standard GPT-style decoder generates operator/variable token sequences conforming to grammar (Valipour et al., 2021).
• Cross-entropy and reward learning: During training, models maximize log-likelihood over token sequences; at inference, top-$k$ sampling or beam search yields candidate expressions, which are post-processed via numerical optimization to fit constants (Biggio et al., 2021).
• End-to-end generative frameworks (e.g., DGSR): Pre-train with loss functions that directly measure normalized mean squared error on data, not just symbolic consistency, thereby collapsing equivalent forms (e.g., accounting for commutativity) (Holt et al., 2023). Inference is accelerated by seeding genetic programming populations with sampler outputs and refining via priority queue training.
• Domain-Aware Symbolic Priors: Incorporate priors derived from mathematical expression corpora specific to physics, biology, or chemistry, integrated as soft regularization (KL divergence) in tree-structured RNN models (Huang et al., 12 Mar 2025).

These models exploit permutation invariance, equation symmetries, and context (via formal grammars). Pre-trained LMs such as SymbolicGPT or DGSR show superior recovery rates, data efficiency, and robustness to variable count compared to both RL-based and classic GP approaches (Valipour et al., 2021, Holt et al., 2023).

4. Search Space Constraints, Complexity Control, and Scalability

Exhaustive search approaches (ESR) and systematic combinatorial enumeration prune the candidate space using:

• Grammar restriction (context-free, e.g., polynomial-of-terms; no nested nonlinearities; controlled maximum depth or variable references) (Kammerer et al., 2021, Desmond, 17 Jul 2025).
• Semantic equivalence and deduplication: Hashing of expression trees modulo algebraic symmetries, folding associativity and commutativity, and constant folding ensure each unique function is considered once.
• Complexity metrics: Node or token counts, description length, or parameter count; exhaustive methods can enumerate all unique expressions up to some maximum complexity $k$.
• Minimum Description Length (MDL): Model selection by joint information-theoretic coding cost: $L_\mathrm{total}(M) = L(M) + L(D|M)$, where $L(M)$ encodes the tree and parameters and $L(D|M)$ is the data misfit (neg log-likelihood) (Desmond, 17 Jul 2025). MDL picks the optimal Pareto point with no arbitrary heuristics.

Systematic search recovers high rates of ground-truth expressions for low complexity and retains partial solutions under noise. However, the exponential scaling in $k$ or number of variables means practical limits must be imposed even for high-performance systems (Kahlmeyer et al., 24 Jun 2025, Kammerer et al., 2021).

5. Extensions: Domain Priors, Experiment Design, and Hybrid Approaches

Incorporating domain knowledge into SR increases tractability and solution relevance:

• Control Variable Genetic Programming (CVGP): Inspired by experimental design, CVGP sequentially “unfreezes” variables under control, fitting reduced forms, then incrementally integrates new degrees of freedom, “freezing” correct subtrees at each step. This yields provable exponential reduction in the search space for certain target families (Jiang et al., 2023).
• Physics-Informed Constraints: Dimensional analysis, invariance to transformations, boundary conditions, or partial known models, can be imposed at the grammar or optimization stage (Austel et al., 2020, Bartlett et al., 17 Dec 2025).
• Neuro-evolutionary and ensemble models: Some approaches evolve symbolic representations of neural networks, with mutation/crossover over topology and gradient-based fine-tuning of parameters. Population diversity and memory-based weight reuse counteract premature convergence (Kubalík et al., 23 Apr 2025).
• Parametric and functional regression: Modern methods extend SR to learning equations with parameter-varying coefficients (parametric PDEs, “hypernetworks”), or treat $g(y) = f(x)$ for jointly searching output transformations and input structure (Zhang et al., 2022, Tohme et al., 2022).

Integration with feature selection and representation learning allows SR to be grafted onto high-dimensional or deep-learning setups as an interpretable “bottleneck” for scientific discovery.

6. Practical Applications and Benchmarking

SR is now routinely applied to:

• Physics: Rediscovery of canonical laws (e.g., Kepler’s, Lorenz, conservation equations), surrogate models in cosmology and BSM physics (AbdusSalam et al., 2024, Bartlett et al., 17 Dec 2025).
• Materials science: Recovery of kinetic laws, free energy functionals (Wang et al., 2019).
• Biology/Ecology: Inference of dynamical systems (Lotka–Volterra, epidemic models) (Brum et al., 27 Aug 2025).
• Engineering: Empirical law identification, surrogate modeling, and control laws.

Benchmarking is standardized on datasets such as the Feynman suite, Nguyen polynomials, and ODE/PDE repositories; core metrics include normalized/test MSE, recovery rate (fraction of runs with exact or equivalent formula), and model complexity (Makke et al., 2022). State-of-the-art methods (Transformer, DGSR, CVGP, domain-prior RNNs) show significant gains in data efficiency, generalizability to higher input dimensions, and noise robustness relative to traditional GP or brute-force search (Holt et al., 2023, Huang et al., 12 Mar 2025, Jiang et al., 2023).

7. Open Challenges and Future Directions

Symbolic regression remains computationally intractable in the worst case. Open areas and trends include:

• High-dimensionality: Scaling SR to larger variable sets and higher complexity relies on hybrid heuristic-genetic, deep generative, or hierarchy-exploiting strategies (Holt et al., 2023, Kahlmeyer et al., 24 Jun 2025).
• Domain-guided grammars: Encoding unit consistency, symmetry, and structure via formal grammars or symbolic priors yields both efficiency and interpretability (Huang et al., 12 Mar 2025, Bartlett et al., 17 Dec 2025).
• Few-shot and zero-shot generalization: Transformer and RNN-based models pre-trained on large synthetic corpora show promise for learning new laws from limited data or extrapolating to unseen regimes (Valipour et al., 2021, Biggio et al., 2021).
• Integration with scientific pipelines: Embedding closed-form surrogates into simulation, optimization, or inference workflows—particularly for "beyond the standard model" physics or material design—delivers orders-of-magnitude speedups (AbdusSalam et al., 2024).
• Model selection and uncertainty: MDL or Bayesian scoring, as in ESR frameworks, provides principled means to resolve the error-complexity trade-off, with applications in automatic model cataloguing (Desmond, 17 Jul 2025).
• Interpretability and reproducibility: Deterministic exhaustive or symmetry-conscious SR enables fully reproducible, white-box model libraries for industrial adoption (Kammerer et al., 2021).

The field continues to advance rapidly in computational efficiency, search-space structuring, and integration of scientific knowledge, making symbolic regression a central component of automated model discovery and scientific machine learning.