Symbolic Extraction of Non-Perturbative Transverse-Momentum-Dependent Distributions from Drell-Yan Data

We present an analytical parametrization of the non-perturbative transverse-momentum-dependent (TMD) parton distribution function of unpolarized quarks, extracted from Drell-Yan data using a combination of neural-network fitting and symbolic regression. A factorized neural network is trained directly against experimental cross-section data from fixed-target, Tevatron, RHIC, and LHC experiments at next-to-next-to-next-to-leading logarithmic accuracy, and symbolic regression is subsequently applied to each network component to discover compact analytical expressions. The final formula is selected from a Pareto front in the space of expression complexity and experimental $\chi^2$, yielding a closed-form non-perturbative function with 9 free numerical constants that achieves $\chi^2/\mathrm{ndf}=1.040$ over 482 data points. A non-trivial $x$-$b_T$ cross term is retained even under a sparsity prior that biases it toward zero, indicating a mild but genuine correlation between the longitudinal momentum fraction and the transverse momentum. This work demonstrates that symbolic regression is a viable tool for bridging flexible machine-learning fits and interpretable analytical TMD parametrizations, and opens a systematic path toward data-driven discovery of specific features of non-perturbative QCD.