Source-linked AI summary

LLM-SR: Scientific Equation Discovery via Programming with Large Language Models

Parshin Shojaee, Kazem Meidani, Shashank Gupta, Amir Barati Farimani, Chandan K Reddy

arXiv:2404.18400v3cs.LGcs.AIcs.CLcs.NE

TL;DR

Equation discovery requires searching a vast space, while symbolic-regression methods often neglect scientific prior knowledge and use constrained representations. LLM-SR treats equations as programs, combining LLM-generated hypotheses with data-driven optimization and evolutionary search. Across custom scientific benchmarks, it outperforms symbolic-regression baselines, particularly on out-of-domain tests.

  • Problem

    Equation discovery seeks accurate, interpretable equations that generalize from observations, but the search space is vast and existing methods often omit scientific prior knowledge or use constrained representations.

    Source: §1 · pp. 1–2
  • Method

    LLM-SR uses an LLM to propose equation programs informed by scientific descriptions, optimizes their parameters against data, and iteratively refines hypotheses through evolutionary search.

    Source: §1 · p. 2, §2.1 · p. 3
  • Results

    LLM-SR consistently outperforms state-of-the-art symbolic-regression methods, discovering more accurate equations with better out-of-domain generalization across its benchmarks.

    Source: §1 · p. 3
  • Takeaways & Limitations

    The authors report that scientific priors make equation-space exploration more efficient, requiring fewer iterations to find accurate equations.

    Source: §1 · p. 3

    Performance depends on the quality and breadth of LLM training data, and iterative LLM queries and parameter optimization may be costly at large scale.

    Source: §5 · p. 10

Abstract

from arXiv · show

Mathematical equations have been unreasonably effective in describing complex natural phenomena across various scientific disciplines. However, discovering such insightful equations from data presents significant challenges due to the necessity of navigating extremely large combinatorial hypothesis spaces. Current methods of equation discovery, commonly known as symbolic regression techniques, largely focus on extracting equations from data alone, often neglecting the domain-specific prior knowledge that scientists typically depend on. They also employ limited representations such as expression trees, constraining the search space and expressiveness of equations. To bridge this gap, we introduce LLM-SR, a novel approach that leverages the extensive scientific knowledge and robust code generation capabilities of Large Language Models (LLMs) to discover scientific equations from data. Specifically, LLM-SR treats equations as programs with mathematical operators and combines LLMs' scientific priors with evolutionary search over equation programs. The LLM iteratively proposes new equation skeleton hypotheses, drawing from its domain knowledge, which are then optimized against data to estimate parameters. We evaluate LLM-SR on four benchmark problems across diverse scientific domains (e.g., physics, biology), which we carefully designed to simulate the discovery process and prevent LLM recitation. Our results demonstrate that LLM-SR discovers physically accurate equations that significantly outperform state-of-the-art symbolic regression baselines, particularly in out-of-domain test settings. We also show that LLM-SR's incorporation of scientific priors enables more efficient equation space exploration than the baselines. Code and data are available: https://github.com/deep-symbolic-mathematics/LLM-SR

1 INTRODUCTION

Symbolic regression seeks interpretable equations that fit observations and generalize, but traditional methods struggle with the search space and often omit scientific prior knowledge. LLM-SR combines LLM-generated equation hypotheses, data-driven optimization, and evolutionary search to discover equations.

1 INTRODUCTION

  • Traditional symbolic-regression methods struggle with the vast equation search space and often omit prior scientific knowledge, motivating methods that integrate such knowledge into the search.Equation discovery seeks models that predict observations, remain interpretable, and generalize to unseen data from the same phenomena.Source: §1 · pp. 1–2
  • LLM-SR iteratively prompts an LLM to propose equation hypotheses, fits them to observed data with optimizers, and uses a dynamic memory of prior equations to guide refinement.It represents equations as Python programs with placeholder parameters, enabling flexible mathematical relations and direct parameter optimization.Source: §1 · pp. 2–3, Figure 1 · p. 2
  • LLM-SR consistently outperforms state-of-the-art symbolic regression methods, finding physically accurate equations with better fit and generalization, especially out of domain.The method explores the equation space more efficiently and requires fewer iterations to find accurate equations.Source: §1 · p. 3
  • The authors evaluate LLM-SR on four custom benchmarks across physics, biology, and materials science, designed to simulate discovery and reduce LLM recitation risk.The paper also reports that natural-language problem descriptions, program representation, data-driven feedback, and iterative refinement are essential components.Source: §1 · p. 3

2 LLM-SR METHODOLOGY

LLM-SR represents equations as programs and uses an LLM to propose equation structures that are optimized against data. An experience buffer of high-scoring programs supplies examples for later iterations.

2.1 PROBLEM FORMULATION

  • Equation program skeletons are represented as functions mapping inputs to targets, with placeholder parameters optimized to maximize data-based fit.The program representation is more expressive than constrained forms such as expression trees, but expands the search space.Source: §2.1 · p. 3

2.2 HYPOTHESIS GENERATION

  • The LLM proposes diverse equation program skeletons from prompts containing instructions, problem specifications, evaluation criteria, and experience examples.Stochastic temperature-based sampling balances exploration and exploitation; programs that fail to execute or exceed a time limit are discarded.Source: §2.2 · p. 4, Figure 2 · p. 4

2.3 HYPOTHESIS OPTIMIZATION AND ASSESSMENT

  • Each skeleton’s placeholder parameters are optimized against observed data, then its fitness is scored as negative mean squared error between predictions and targets.The method uses numpy+BFGS for problems with fewer parameters and torch+Adam for larger-scale problems benefiting from differentiable programming.Source: §2.3 · pp. 4–5

2.4 EXPERIENCE MANAGEMENT

  • Experience management stores diverse high-quality programs and samples them as in-context examples to guide subsequent LLM proposals.The framework iterates through example sampling, prompt updating, skeleton generation, and evaluation, adding hypotheses when they improve on the best score.Source: §2.4 · p. 5

3 EXPERIMENTS

The experiments evaluate LLM-SR on custom scientific benchmarks and compare it with symbolic-regression baselines. They report stronger out-of-domain performance, more interpretable equations, and component ablations that identify contributions from prior knowledge and iterative optimization.

3.1 BENCHMARKS AND DATASETS

  • The four custom benchmarks span physics, biology, and materials science, using modified models, nonlinear relations, and experimental data to limit recitation and simulate discovery.The oscillator designs use nonlinear structures beyond common systems; the stress-strain task has no predetermined theoretical model structure.Source: §3.1 · p. 6

3.3 QUANTITATIVE RESULTS

  • On E. coli growth, LLM-SR reaches OOD NMSE ∼0.0037 versus > 1 for other methods, while outperforming baselines across benchmarks with 2.5K versus 2M+ iterations.Table 1 reports NMSE on in-domain and out-of-domain test sets; the performance gap is more pronounced out of domain.Source: Table 1 · p. 7, §3.3 · p. 7, Figure 3 · p. 7
  • LLM-SR variants navigate the equation search space more quickly, while symbolic-regression baselines fail to match their performance even after 2M+ iterations.The authors attribute the faster search to scientific prior knowledge and LLMs acting as mutation or crossover agents in iterative refinement.Source: §3.3 · pp. 7–8, Figure 3 · p. 7

3.4 QUALITATIVE ANALYSIS

  • For the oscillator problems, LLM-SR recovers more true symbolic terms and combines interpretable driving, damping, and restoring-force terms; on E. coli growth, its predicted distribution aligns with ground truth in OOD regions.PySR and uDSR deviate significantly from ground truth in OOD regions, while the baselines’ oscillator equations lack clear physical interpretation.Source: §3.4 · p. 8, Figure 4 · p. 8, Figure 5 · p. 9

3.5 ABLATION STUDY

  • On Oscillation 2, removing scientific prior knowledge or iterative refinement substantially harms performance, and generating full equations without separate parameter optimization also worsens results.Removing program representation has a negative but less severe effect; numpy+BFGS performs slightly better than torch+Adam, which the authors attribute most likely to LLM code-generation proficiency.Source: §3.5 · p. 9, Figure 6 · p. 9

4 RELATED WORK

Prior symbolic-regression methods span search-based, learning-based, and hybrid approaches, but often lack scientific priors or expressive equation representations. LLM-SR applies LLM scientific knowledge and code generation to equation discovery with data-driven evaluation.

LLMs and Optimization

  • LLM-SR combines LLM scientific priors and code generation with data-driven evaluators, placing equation discovery within iterative optimization approaches that use feedback.The related work identifies FunSearch as a system combining LLMs with systematic evaluators to search for mathematical heuristics.Source: §4 · pp. 9–10

LLMs for Scientific Discovery

  • Despite growing LLM use in scientific tasks, equation discovery and symbolic regression remain relatively unexplored applications of LLM scientific knowledge and reasoning.Prior work has used LLMs to propose hypotheses and, with tools or simulators, support statistical discovery and reasoning.Source: §4 · p. 10

Symbolic Regression

  • Search-based symbolic regression explores equation structures and parameters but can struggle with scalability and efficiency; learning-based models map observations to expressions using synthetic data and Transformers.Hybrid methods use neural priors to guide search, aiming to combine expressiveness and efficiency.Source: §4 · p. 10

5 CONCLUSION AND FUTURE WORK

LLM-SR uses LLM scientific knowledge and code generation with evolutionary search to discover equations, with reported gains particularly on out-of-domain tests. The authors note limits from LLM training data and computational cost; the appendix describes the evaluated SR baselines.

5 CONCLUSION AND FUTURE WORK

  • LLM-SR treats equations as programs and combines LLM-generated hypotheses with evolutionary search, achieving superior benchmark performance particularly in out-of-domain settings.The authors identify training-data breadth and the cost of iterative LLM queries and parameter optimization as limitations.Source: §5 · p. 10

A.1 MODELS

  • The baseline set spans genetic programming, reinforcement-learning symbolic regression, hybrid search, and pretrained Transformer approaches.The listed methods include GPlearn, PySR, DSR, uDSR, NeSymReS, and E2E.Source: §A.1 · p. 14

A.2 DATA PREPROCESSING AND MODEL EXECUTION

The appendix describes baseline execution and illustrates prompts for two benchmark problems. Transformer SR models receive normalized inputs, while search-based baselines run for more than two million iterations.

A.2 DATA PREPROCESSING AND MODEL EXECUTION

  • Transformer SR inputs are normalized; search-based baselines run for over 2 million iterations, with five replications and the best results reported.The stated baselines are GPlearn, PySR, DSR, and uDSR.Source: §A.2 · p. 14
  • The prompt examples show E. coli Growth and Stress-Strain problem specifications with a simple linear equation skeleton; evaluation-function details are omitted.The caption directs readers to Fig. 2 for those details.Source: Figure 7 · p. 15

B DETAILS OF LLM-SR METHOD AND IMPLEMENTATION

LLM-SR prompts an LLM with problem context and accumulated equation examples, evaluates proposed programs against data, and uses an island-based experience buffer to guide later generations.

Hypothesis Generation and Data-driven Evaluation

  • Natural-language problem context and differentiable PyTorch specifications guide the LLM in generating scientifically plausible equation programs.The prompt encourages smooth approximations for non-differentiable operations such as if-else conditions.Source: B DETAILS OF LLM-SR METHOD AND… · p. 15, Figure 8 · p. 16

Implementation

  • Each prompting step generates four equation skeletons, with at most ten parameters, evaluated concurrently under 30-second and 2-GB limits.The setup uses four evaluators and discards programs that exceed either resource limit.Source: B DETAILS OF LLM-SR METHOD AND… · pp. 15–16

Experience Buffer Management

  • Ten islands store candidate programs; periodic resets replace the worst-performing half with copies of high-scoring programs from surviving islands.Programs are clustered by score signature to retain candidates with different performance characteristics.Source: B DETAILS OF LLM-SR METHOD AND… · pp. 16–17

Experience Sampling

  • Prompts sample examples by first selecting an island, then favoring higher-scoring clusters and shorter programs within the chosen cluster.The sampled programs become in-context demonstrations for generating new candidates.Source: B DETAILS OF LLM-SR METHOD AND… · p. 17

C LIMITATION OF FEYNMAN BENCHMARK PROBLEMS

The authors argue that LLMs may recall familiar Feynman equations rather than discover them, making this widely used benchmark an imperfect test of equation discovery.

Benchmark limitation

  • The 120 fundamental-physics problems are widely used, but likely memorization may keep them from reflecting LLMs’ true discovery capabilities.The perplexity analysis measures equation likelihood conditional on scientific context.Source: C LIMITATION OF FEYNMAN BENCHMARK… · p. 17

Perplexity Analysis

  • Feynman problems have significantly lower perplexity than the new benchmarks, indicating greater prediction certainty and a higher likelihood of recitation.Perplexity is computed for equation generation given scientific context; Stress-Strain is excluded because its experimental data lack a predetermined theoretical model structure.Source: Figure 9 · p. 17, C LIMITATION OF FEYNMAN BENCHMARK… · p. 18

Discovery Error Curve Analysis

  • LLM-SR reaches low NMSE on Feynman problems within very few iterations, often one pass, while the new benchmarks require reasoning and exploration.One-pass responses often recite the exact corresponding physics expressions; the trajectory compares best NMSE against iterations.Source: Figure 11 · p. 18, C LIMITATION OF FEYNMAN BENCHMARK… · p. 18

D ADDITIONAL DETAILS ON NEW BENCHMARK PROBLEMS

The section describes custom benchmarks designed to challenge equation discovery with novel oscillator and bacterial-growth models, alongside generated datasets for those problems.

D ADDITIONAL DETAILS ON NEW BENCHMARK PROBLEMS

  • Oscillation 1, Oscillation 2, and E. coli Growth data were generated for this work, while the Aluminium 6061-T651 stress-strain dataset is publicly available.The generated datasets are released with the LLM-SR repository, and the stress-strain data carry a CC BY 4.0 license.Source: D ADDITIONAL DETAILS ON NEW BENCHMARK… · p. 18

D.1 NONLINEAR OSCILLATOR EQUATIONS

  • Two custom nonlinear oscillators extend beyond common oscillator systems with diverse nonlinear structures, creating challenging but solvable and physically plausible discovery tasks.Their simulated trajectories are partitioned into training, in-domain validation, and out-of-domain validation, with T = [0, 20) used for out-of-domain evaluation.Source: §D.1 · pp. 19–20

D.2 E. COLI GROWTH RATE EQUATIONS

  • The E. coli benchmark combines a linear population-density term and Monod substrate prior with novel temperature and pH dependencies.The temperature function combines hyperbolic tangent growth acceleration with quartic high-temperature decay; the custom pH function combines exponential and trigonometric terms.Source: §D.2 · p. 21
  • The custom temperature and pH functions preserve characteristic bacterial-growth responses while narrowing optimal ranges and steepening declines relative to established models.The custom benchmark combines biological prior knowledge with novel function/operator combinations to challenge equation discovery without making the relationships biologically implausible.Source: §D.2 · p. 21, Figure 15 · p. 22

D.3 MATERIAL STRESS BEHAVIOR ANALYSIS

The section presents Aluminium 6061-T651 stress-strain data as an empirical equation-discovery case without a universally accepted theoretical model, and reports LLM-SR’s comparative performance across benchmark tasks.

D.3 MATERIAL STRESS BEHAVIOR ANALYSIS

  • The experimental dataset covers uniaxial tensile behavior of Aluminium 6061-T651 at six temperatures from 20°C to 300°C, with 200°C held out for out-of-domain validation.The curves show decreasing yield strength and ultimate tensile strength as temperature rises, and generally increasing strain at failure.Source: §D.3 · p. 22
  • Because stress-strain-temperature relationships for specific materials often lack a universally accepted theoretical model, the problem calls for exploratory, data-driven empirical equation discovery.The study uses this real-world case to assess discovery of physically meaningful, interpretable equations where empirical modeling often dominates.Source: §D.3 · pp. 22–23

Ablation comparison

  • LLM-SR consistently outperforms PySR and its ablated variants on in-domain and out-of-domain errors across four benchmarks.The comparison covers Oscillation 1, Oscillation 2, E. coli Growth, and Stress-Strain.Source: Figure 17 · p. 23

E ADDITIONAL ABLATION STUDY

Ablations show that both scientific prior knowledge and program representation contribute to LLM-SR’s performance. The multi-island buffer and sampling design also improves performance in both test settings by supporting exploration.

E ADDITIONAL ABLATION STUDY

  • Across all benchmarks and both test settings, removing either prior knowledge or program representation underperforms the full LLM-SR model.The two components have comparable effects on both oscillation problems; prior knowledge has a stronger effect for E. coli growth and stress-strain OOD performance.Source: E ADDITIONAL ABLATION STUDY · p. 23
  • For the oscillation problems, prior knowledge and code generation have comparable impacts on model performance in both test settings.For E. coli growth, prior knowledge has a more pronounced effect; for stress-strain, its effect is stronger on OOD performance.Source: E ADDITIONAL ABLATION STUDY · p. 23
  • Fewer buffer islands reduce exploration, diversify hypotheses less, and can cause premature convergence to early equation structures.The single-island variant uses deterministic top-k selection for in-context examples; the multi-island design improves performance in both in-domain and OOD settings.Source: E ADDITIONAL ABLATION STUDY · p. 24, Table 3 · p. 24

F.1 NOISE ROBUSTNESS ANALYSIS

The supplementary analyses examine noise robustness, comparisons with LLM-based optimization methods, and a larger LLM backbone. LLM-SR is more resilient to noise, outperforms LMX and FunSearch, and GPT-4o yields slightly better performance than GPT-3.5 and Mixtral on the tested oscillator datasets.

F.1 NOISE ROBUSTNESS ANALYSIS

  • At moderate Gaussian noise levels (σ = 0.01, σ = 0.05), LLM-SR maintains lower NMSE than PySR, particularly for OOD predictions.Both methods degrade as noise increases on Oscillation 2, but LLM-SR shows better noise resilience.Source: Figure 18 · p. 24, §F.1 · p. 24

F.2 LLM-BASED OPTIMIZATION BASELINES

  • LLM-SR and its variants outperform LMX and FunSearch in both in-domain and OOD settings on the oscillator problems.The full LLM-SR achieves the best overall performance; traditional SR baselines also outperform the existing LLM-based approaches on both datasets.Source: §F.2 · p. 25, Table 4 · p. 25

F.3 LLM BACKBONES

  • LLM-SR with GPT-4o performs slightly better than with GPT-3.5 and Mixtral on Oscillation 1 and 2, particularly in OOD settings.The GPT-4o comparison used the same 2500 iterations; computational budget limited these experiments to the two oscillation datasets.Source: §F.3 · p. 25, Table 5 · p. 25

G ADDITIONAL QUALITATIVE RESULTS

The qualitative analyses trace how LLM-SR’s equation programs evolve across oscillation, growth, and stress-strain problems. They describe domain-informed search choices, an approximate match to the Oscillation 2 equation over the observed displacement range, and prediction alignment with ground truth.

Oscillation 2

  • For Oscillation 2, LLM-SR finds common nonlinear terms early and later discovers a form whose 5(1 − exp(x)) term approximates −5x·exp(0.5x) for x in (−2, 2).That displacement interval is the approximate range in the dataset.Source: G ADDITIONAL QUALITATIVE RESULTS · p. 26

Oscillation 1

  • On Oscillation 1, the evolving skeletons incorporate candidate nonlinear terms associated with driving, restoring, and damping forces.These terms appear in code comments or variable names as the model aims to improve accuracy.Source: G ADDITIONAL QUALITATIVE RESULTS · p. 26, Figure 20 · p. 27

E. coli growth

  • For E. coli growth, LLM-SR explores temperature and pH effects, recognizing potential optima early and combining function components in the final step.The account attributes recognition of bell-shaped temperature and pH effects to model prior knowledge; search focuses more on these variables than substrate concentration.Source: G ADDITIONAL QUALITATIVE RESULTS · p. 26, Figure 21 · p. 28

Stress-Strain

Discovered Equations

Behavior of Discovered Models

Loading 2404.18400v3…