PaperlayerView exact revision on arXiv

Source-linked AI summary

Symbolic Density Estimators for Unnormalized Distributions

Vikas Kanaujia, Riyansha Singh, Shashank Sharma, Vipul Arora

arXiv:2610.10807v1cs.LGhep-lat

TL;DR

Inferring symbolic probability densities from data is challenging because analytical forms have traditionally depended on domain expertise. The paper combines deep generative models with symbolic regression and physical inductive biases; experiments support estimating Hamiltonian expressions for toy distributions and lattice systems. The framework also produces compact symbolic approximations for renormalized ϕ4 theory.

  • Problem

    Inferring analytical or symbolic probability-density forms directly from samples is challenging and has traditionally relied on expert-selected forms and parameter estimation.

    Source: §1 · pp. 1–2

  • Method

    The framework uses likelihood- or score-based deep generative models to supply targets for symbolic regression, with factorization, task-specific constraints, and noisy perturbations.

    Source: §1 · pp. 2–3

  • Results

    Experiments support estimating Hamiltonian expressions for the XY model and ϕ4 theory, density functions for toy distributions, and compact symbolic approximations for renormalized ϕ4 theory.

    Source: §6 · p. 16

  • Takeaways & Limitations

    The framework provides interpretable symbolic equations for locally interacting systems while using domain-specific prior knowledge.

    Source: §7 · p. 16

    Expression quality depends on the chosen symbolic library and patch decomposition, so the framework is not a fully automated discovery system.

    Source: §5 · p. 15

Abstract

from arXiv · show

Estimating the symbolic or analytical form of probability density functions (PDFs) from observed samples is a fundamental challenge in statistical and computational modelling. This process is critical for deriving interpretable and generalizable relationships characterizing the underlying phenomenon. Traditionally, this estimation depends strongly on domain expertise and prior field-specific knowledge, with experts selecting appropriate functional forms or parametric families based on empirical evidence and theoretical understanding. The coefficients of these forms are then typically determined through parameter estimation. In this paper, we develop a framework for estimating symbolic expressions of unnormalized distributions from observed samples using domain-specific prior knowledge, such as the range of interactions and a predefined set of primitive functions. We integrate deep generative models with symbolic regression (SR), incorporating inductive biases, such as factorizing large distributions, to keep the problem tractable. The deep generative models we examine include likelihood-based models, viz., flow models, and score-based models. Experiments show the effectiveness of the proposed framework for estimating density functions for multivariate toy distributions as well as lattices from computational physics, namely, XY model and $\phi^4$ theory. When applied to the renormalization problem in $\phi^4$ theory, the proposed framework estimates compact symbolic approximations of the hamiltonian function at different scales directly from samples, yielding expressions that may be challenging to derive using traditional perturbative or analytic approaches in nonperturbative settings.

1 Introduction

The paper addresses the challenge of inferring symbolic density representations without relying on experts to choose their analytical form. It combines deep generative models with symbolic regression and physical inductive biases, applying the framework to toy distributions and lattice systems.

1 Introduction

  • The framework targets the labor-intensive, expert-dependent practice of selecting analytical probability models and estimating their parameters from data.Automatically inferring a symbolic density representation directly from observations is described as a significant challenge.Source: §1 · pp. 1–2
  • Existing deep generative approaches provide likelihoods or scores but are black-box estimators, while many symbolic regression methods work well only with few variables.Prior attempts to recover low-dimensional unnormalized distributions rely on assuming the distribution’s form.Source: §1 · p. 2
  • Applications include multivariate Gaussian and many-well distributions, the XY model, scalar ϕ4 theory, and renormalized ϕ4 effective Hamiltonians.The paper highlights the renormalization setting as analytically challenging in non-perturbative conditions.Source: §1 · p. 2
  • The framework integrates deep generative models with symbolic regression to estimate interpretable energy functions or Hamiltonians from samples of unnormalized distributions.Likelihood or score information from learned models converts density estimation into a supervised symbolic regression problem.Source: §1 · p. 2
  • Local Hamiltonian factorization, task-specific symbolic constraints, and multi-scale noise perturbations address dimensionality, prior knowledge, and low-density sample coverage.Factorization exploits spatial locality; perturbations augment high-density samples with controlled noise to improve low-density coverage.Source: §1 · pp. 2–3

2 The Framework

The framework first trains a deep density or score model, then uses its likelihood or gradient information as targets for symbolic regression. It incorporates local factorization, perturbations, and alternative score-based formulations to make symbolic estimation more tractable.

2 The Framework

  • Training proceeds from a density or score model to symbolic regression on samples and model-provided likelihood or gradient information, yielding an energy expression.Density factorization can conform to known system structure.Source: §2 · pp. 3, 6
  • Symbolic regression searches a library of primitive operators for an explicit expression that fits input-output data while balancing accuracy and interpretability.Unlike conventional regression with a predefined model, it searches for both equation structure and parameters.Source: §2 · pp. 3–4
  • A trained flow supplies negative log-likelihood targets for symbolic regression, and terms independent of x are dropped to estimate the energy.The flow is trained from target-distribution samples by minimizing KL divergence.Source: §2 · p. 4
  • Score matching estimates the gradient of log density, which is independent of the normalizing constant, and symbolic regression can recover energy by integrating the fitted field.This recovery requires a conservative field; approximation errors can produce small deviations in practice.Source: §2 · p. 5
  • Denoising score matching perturbs data across noise levels to improve score estimation in low-density regions, while an integration-free variant models the energy directly.The direct-energy formulation obtains the score as the negative gradient of a neural energy model, avoiding symbolic integration of a fitted score expression.Source: §2 · p. 5
  • Local lattice factorization represents total energy through site-wise local energies, and patch-based symbolic regression estimates each local function over fewer variables.The patch structure can reflect interactions with any number of neighbors.Source: §2 · pp. 5–6
  • Patch-based energy representations can be redundant, whereas fitting site-wise gradients avoids that representation redundancy.Different local functions can yield the same total energy under the factorized representation.Source: §2 · p. 6
  • Noisy perturbations augment observed samples from low-energy regions so symbolic regression also receives targets from high-energy regions.The deep models provide likelihood or score estimates for the perturbed samples.Source: §2 · p. 6

3 Case Studies

The case studies apply the framework to toy distributions and lattice systems, then estimate effective Hamiltonians from coarse-grained samples in non-perturbative ϕ4 theory.

Framework

  • The estimation pipeline trains a generative model, augments low-density samples, constructs labels, and uses symbolic regression to recover an analytical expression.For locally factorisable distributions, symbolic regression operates on local patches.Source: §3 · p. 7

Multivariate Gaussian

  • For the multivariate Gaussian case, the study uses two-dimensional distributions with a fixed mean and diagonal or full covariance variants, expecting recovery of their Hamiltonians.The two variants differ in whether the covariance is diagonal or full.Source: §3 · p. 7

Many Well Distribution

  • The many-well study evaluates Hamiltonian estimation across dimensions d = 2, 4, 8, 16, 32, and 64.The many-well distribution is constructed from d/2 copies of a two-dimensional double-well distribution.Source: §3 · pp. 7–8

Lattice Systems

  • The lattice case studies include the XY model on an 8 × 8 grid with λ = 1/1.4 and scalar ϕ4 theory on 8 × 8 lattices.A 2 × 2 patch is used to factorize the XY Hamiltonian; for ϕ4, a 2 × 2 patch suffices, while larger patches capture distant interactions but make symbolic regression harder.Source: §3 · p. 8

Renormalization

  • The framework estimates effective Hamiltonians directly from coarse-grained samples, where analytic methods are limited and perturbative analysis fails.For ϕ4 theory with λ1 = 4 and λ2 = 1, samples are successively coarsened from 64 × 64 to 32 × 32, 16 × 16, and 8 × 8 lattices.Source: §3 · p. 8

4 Experiments and Results

Experiments test symbolic density estimation across toy distributions, lattice models, and renormalized ϕ4 systems. Results include agreement with coarse-grained observables and recovery of local interaction structure, with performance varying by model and symbolic-regression method.

Evaluation

  • Evaluation separates predictive accuracy, assessed with MSE and R2, from symbolic-expression quality because similar predictive scores can accompany structurally different expressions.The symbolic assessment also considers coefficient accuracy, structural correctness, and spurious interactions.Source: §4 · pp. 9–10

Multivariate Gaussian

  • For both diagonal and non-diagonal 2-D Gaussians, the methods recover the true H(x); PySR achieves performance comparable to EQL using 10^3 rather than 2.5 × 10^4 training samples.MeSSY achieves MSE close to zero and R2 near one but produces incorrect coefficients and additional spurious terms.Source: §4 · p. 10, Table 1 · p. 10

Many-Well distribution

  • For Many-Well distributions, both model families and SR methods recover the true H(x) at d = 2, 4, and 8, while denoising score matching beats flows at d = 16, 32, and 64.At higher dimensions, PySR also surpasses EQL; MeSSY does not converge for d > 2.Source: §4 · p. 11, Table 2 · p. 11

XY model

  • For the XY model, the flow model outperforms score-based methods, denoising score matching beats standard score matching, and PySR outperforms EQL; MeSSY fails to converge.Table 3 presents the estimated expressions against the true XY Hamiltonian, with incorrect terms marked in gray.Source: §4 · p. 11, Table 3 · p. 12

Scalar ϕ4 theory

  • For ϕ4 theory, PySR recovers the exact local Hamiltonian across patch sizes without spurious terms, whereas EQL adds 4, 4, and 7 terms for 2 × 2, 3 × 3, and 4 × 4 patches.Larger patches expand the symbolic search space, but PySR retains only physically relevant nearest-neighbor interactions.Source: §4 · p. 12, Table 5 · p. 13

Renormalization

  • Histogram overlap exceeds 98% across observables and lattice scales, while EMD ranges from approximately 10^-4 to 10^-3.The comparison uses magnetization, absolute magnetization, and magnetic susceptibility from samples generated with estimated Hamiltonians versus coarse-grained data.Source: §4 · p. 13, Table 6 · p. 14
  • As lattice size decreases, the estimated x4 coefficient becomes negligible and is dropped, while observable means and standard deviations closely match coarse-grained data.The coefficient trend also appears in the non-perturbative regime; at 32 × 32, mean absolute magnetization is 0.0053 for coarse-grained data and 0.0055 for generated data.Source: §4 · pp. 13–14, Table 7 · p. 14, Figure 3 · p. 14

5 Discussion

The discussion reports lower MSE after symbolic regression (SR) than after the flow model, while identifying generative-model accuracy, sample coverage, SR scalability, and domain choices as constraints. It also notes a score-estimation challenge for the XY model's periodic state space.

5 Discussion

  • The two-stage pipeline learns a generative model first and then applies SR; the authors suggest low-complexity SR may reduce noise in flow-model likelihood information.The authors report lower MSE after SR and describe noise reduction as a plausible explanation.Source: §5 · pp. 14–15
  • Performance depends on generative-model accuracy, state-space coverage, and SR scalability; for highly multimodal targets, flow mode-covering can bias estimates and yield coefficient errors and spurious terms.The passage links these effects to larger errors in estimated Hamiltonians.Source: §5 · p. 15
  • Adequate sample coverage is an assumption; in MW-64, about 232 metastable modes make representative sampling difficult and symbolic-expression accuracy is reduced.The passage notes that under-sampling important state-space regions can degrade both density estimation and symbolic expression.Source: §5 · p. 15
  • SR search grows rapidly with variable count, and long-range interactions require larger patches that expand the search space and can increase coefficient errors and spurious terms.The discussion says the SR component often constrains overall scalability more than generative modeling does.Source: §5 · p. 15
  • The framework relies on domain knowledge for its symbolic library and patch structure, so excluded forms or interactions can leave only approximate Hamiltonian descriptions.The estimated expressions are characterized as compact, interpretable approximations rather than exact representations of the physical system.Source: §5 · p. 15
  • For the XY model, Gaussian perturbations and Euclidean score fields do not fully respect its periodic manifold, potentially increasing score errors near boundaries and accumulating error during integration.The passage identifies oscillatory trigonometric interactions as another challenging region for score estimation.Source: §5 · p. 15
  • MSE is lower after SR with PySR or EQL than after the flow model, according to the reported comparison.Table 8 compares MSE from the flow model with MSE after SR using flow-model outputs.Source: §5 · p. 15, Table 8 · p. 16
  • In challenging settings, the estimated expressions should be read as compact symbolic approximations of effective Hamiltonians, not exact forms of the underlying dynamics.The discussion states that this distinction applies in more challenging settings.Source: §5 · p. 16

6 Conclusion

The conclusion presents a framework combining deep generative models and symbolic regression to estimate Hamiltonians from samples using domain-specific priors. It reports applications to toy distributions and lattice systems, including scale-consistent symbolic approximations in non-perturbative renormalization, while noting limits from SR scalability and prior choices.

6 Conclusion

  • The framework combines deep generative models and SR, using density factorization and sample perturbation to estimate Hamiltonians from samples with domain-specific priors.The conclusion describes likelihood-based and score-based generative models alongside symbolic regression.Source: §6 · p. 16
  • Experiments estimate expressions for toy distributions, XY and ϕ4 systems, and compact effective Hamiltonians consistent with coarse-grained dynamics across scales in non-perturbative renormalization.The conclusion highlights cases where exact coarse-grained Hamiltonians are unknown and analytically challenging in the considered setting.Source: §6 · p. 16
  • Existing SR methods' limits with very large variable counts hinder discovery of long-range or non-local lattice interactions, while expression quality depends on symbolic libraries and patch structures.The conclusion reports good results with patches up to 4 × 4 (16 variables).Source: §6 · p. 16

7 Impact Statement

The impact statement presents interpretable symbolic equations from data as a way for scientists to use machine learning while analyzing large datasets. It points to possible applications across computational sciences and to challenging systems whose equations are difficult to derive analytically.

7 Impact Statement

  • The framework's potential applications include Physics, Chemistry, and Climate science.The statement identifies these as areas in computational sciences.Source: §7 · p. 16
  • The approach aims to let scientists analyze large datasets with machine learning while estimating interpretable symbolic equations for locally interacting systems using domain-specific priors.The statement frames this as addressing concern about machine learning's ability to provide insight into studied phenomena.Source: §7 · p. 16
  • The statement highlights systems with governing equations difficult to derive analytically, for which approximate symbolic forms can still be estimated from observed data.It specifically notes the considered non-perturbative setting.Source: §7 · pp. 16–17
  • Suggested follow-up areas include phase transitions and renormalization, nuclear interaction potentials from scattering data, and climate extremes from weather data.Source: §7 · p. 17

A Data Generation

The appendix describes sample-generation protocols for Gaussian and many-well toy distributions and for XY and scalar ϕ4 lattice systems.

A Data Generation

A.1 Multivariate Gaussian Distribution

  • The Gaussian experiments used 25,000 training, 5,000 validation, and 5,000 testing samples, plus 20,000 separate inference samples; a full-covariance setting is also mentioned.The first setting has mean [−1.5, 1.5] and unit covariance.Source: §A.1 · p. 19

A.2 Many-Well Distribution

  • Many-Well samples were generated by rejection sampling for dimensions 2, 4, 8, 16, 32, and 64, with 100,000 training, 10,000 validation, 10,000 testing, and 20,000 inference samples per variant.Source: §A.2 · p. 19

A.3 XY Model distribution

  • The XY dataset used Metropolis-Hastings sampling on 8×8 lattices, discarded about 30,000 burn-in samples, and retained configurations every 320 steps to reduce training-sample correlation.The passage reports 100,000 training, 10,000 validation, and 10,000 testing samples.Source: §A.3 · p. 20

A.4 Scalar ϕ4 theory distribution

  • Scalar ϕ4 samples used HMC on 8×8 lattices for mass settings λ2 = {−4, 1} with λ1 fixed at 4; the renormalization experiment used 64×64 lattices with λ1 = 4.0 and λ2 = 1.The 8×8 runs discarded the first 1,000 samples and used 100,000 training, 10,000 validation, and 10,000 testing samples per mass setting; the 64×64 run discarded 5,000 and used 100,000 training and 25,000 each for validation and testing.Source: §A.4 · p. 20

B Model and Training details

The experiments pair flow- or score-based density models with symbolic regression, using benchmark-specific architectures and physically motivated search constraints. These choices make the hypothesis space reflect prior knowledge about interactions and Hamiltonian terms.

B.1 Normalizing Flows

  • Normalizing flows use RNVP coupling layers; the XY data are transformed from a circular manifold into Euclidean space to suit the Gaussian base distribution.The flow setup varies by benchmark, including convolutional layers for lattice data.Source: §B.1 · pp. 20–21

B.2 Score Model

  • Score models estimate the score function from data, using score matching for toy distributions and denoising score matching for XY and scalar ϕ4.The score architecture follows Song & Ermon (2019), adapted for the lattice datasets.Source: §B.2 · p. 21

B.3 EQL

  • EQL embeds symbolic operators in neural networks for end-to-end training and uses smoothed L0.5 regularization to encourage sparse, interpretable expressions.Primitive sets vary by benchmark; for ϕ4, quartic and product operators are restricted or excluded to constrain the search.Source: §B.3 · pp. 21–22
  • XY-model EQL adds sine and cosine features to handle circular topology while avoiding nested trigonometric functions that complicate score integration.Inputs are expanded into sine and cosine variants before symbolic estimation.Source: §B.3 · p. 22

B.4 PySR

  • PySR uses a multi-population tree-based genetic algorithm, with sample size and maximum expression complexity selected using training-loss convergence.Complexity is capped because further increases exponentially raise search time.Source: §B.4 · p. 22

C Inductive Biases and Symbolic Search Space.

  • Symbolic search is constrained by operator libraries, expression complexity, and patch sizes chosen to cover known interaction ranges, so estimates are not assumption-free.Quartic, multiplicative, and trigonometric operators encode expected physics; patch sizes match or exceed the known interaction range.Source: C Inductive Biases and Symbolic Search… · pp. 22–23

D Impact of Hyperparameter Choices on Method Performance

The section examines how training data, model complexity, and learning rate affect the generative and symbolic-regression stages, and explains how reported errors are normalized.

D Impact of Hyperparameter Choices on Method Performance

  • PySR is more sensitive to training-sample count and complexity, whereas EQL is more sensitive to learning rate.Generative-model hyperparameters are selected using validation performance; SR hyperparameters are adjusted against validation error and model complexity.Source: D Impact of Hyperparameter Choices on… · p. 23
  • EQL typically needs 10k–100k samples for symbolic regression, versus PySR's 500–2000; requirements rise with dimensionality and mode count.MW-64 requires more samples than Φ4 (d = 64) to achieve low complexity and low MSE because it has more modes.Source: D Impact of Hyperparameter Choices on… · p. 23
  • When model complexity is too high or samples miss distribution modes, both PySR and EQL can produce uninterpretable equations; generative-model errors may propagate into spurious terms.For PySR, more data can initially suppress extra terms, but complexity eventually dominates; EQL also remains sensitive to model complexity and upstream errors.Source: D Impact of Hyperparameter Choices on… · p. 24

E Metrics

  • Because Hamiltonians from flow- and score-based methods are defined up to an additive constant, evaluation removes the mean offset before calculating MSE and R2.The same normalization is applied after integrating a predicted score field, which introduces an arbitrary integration constant.Source: E Metrics · pp. 24–25

F Robustness Analysis

The robustness analysis repeats symbolic-recovery experiments across random seeds for scalar ϕ4 and many-well benchmarks, using per-seed results and aggregate distributions to assess variability.

SM + PySR (18 seeds)

  • Five reported cases had MSE 0.2–0.5 and R2 >0.99, while three had MSE 0.5–1.0 and R2 0.98–0.99.Five cases had MSE 1.0–4.0 and R2 0.90–0.98, and five had MSE >4.0 and R2 <0.90.Source: SM + PySR (18 seeds) · p. 26

F. Robustness Analysis

G Additional Experiments on Long-range Interaction Systems

Additional experiments test long-range interactions and symbolic recovery when the primitive library lacks target terms. Results show successful estimation under some constrained settings, with accuracy depending on patch size and symbolic-regression choices.

G Additional Experiments on Long-range Interaction Systems

  • With sufficiently large patches, NF + PySR and DSM + PySR recover long-range interaction structure with perfect precision and recall and no spurious terms.EQL pipelines retain many spurious terms, particularly NF + EQL, despite relatively high recall.Source: G Additional Experiments on Long-range… · p. 28, Table 16 · p. 28
  • Long-range interactions require patches covering the full interaction neighborhood, but larger patches expand the SR search space and can reduce accuracy.The authors identify symbolic regression, rather than generative modeling, as the primary challenge in these settings.Source: G Additional Experiments on Long-range… · pp. 27–29

H Symbolic Recovery Under Primitive Library Misspecification

  • In the ϕ4 experiment without explicit polynomial primitives, PySR achieves R2 = 0.99 and estimates the dominant interaction structure, while EQL degrades and adds spurious terms.The target remains representable through repeated multiplication, so this tests indirect recovery rather than strict non-representability.Source: H Symbolic Recovery Under Primitive… · p. 29, Table 17 · p. 30

H.1 Stronger Misspecification: XY Model Without Trigonometric Operators

  • Without trigonometric operators, PySR and EQL approximate the XY cosine interaction with polynomial compositions, achieving R2 values of 0.81 and 0.85, respectively.The resulting expressions capture low-order terms but do not exactly match the true Hamiltonian.Source: §H.1 · pp. 30–31, Table 18 · p. 30
  • Under primitive-library restrictions, useful approximations remain possible, but expressions should be treated as effective approximations rather than exact Hamiltonian recovery.The paper links expression quality to both primitive-set expressiveness and the regressor's ability to compose available operators.Source: §H.1 · pp. 30–31

H.2 Error Propagation Under Primitive Set Misspecification

The section tests whether symbolic regression reduces approximation error in learned generative models when the primitive library is restricted. It finds that error reduction depends on how well the symbolic regressor can represent the target within the available hypothesis space.

H.2 Error Propagation Under Primitive Set Misspecification

  • The analysis compares learned-flow and symbolic-Hamiltonian errors against the true Hamiltonian, reported in Table 19.MSE(H, Hθ) measures learned-flow error; MSE(H, ˆH) measures symbolic-approximation error.Source: §H.2 · p. 31, Table 19 · p. 31
  • 0.22 to 0.01: in the ϕ4 experiment, PySR reduces MSE relative to the learned flow model, while EQL adds spurious terms and yields larger error.PySR estimates the missing polynomial structure through compositions of the remaining primitives.Source: §H.2 · p. 31
  • In the XY experiment, both PySR and EQL reduce error, with EQL achieving lower symbolic estimation error, but the reduction is smaller than in ϕ4.Removing trigonometric operators leaves the cosine interaction unrepresentable directly, so it must be approximated through polynomial compositions.Source: §H.2 · p. 31
  • Error reduction depends on primitive-set misspecification and SR expressiveness: it can be substantial when the target remains representable, but is more limited when its functional form is absent.Under stronger misspecification, symbolic regression can still produce meaningful approximations.Source: §H.2 · p. 32
Loading 2610.10807v1…