Reconstruction of dark energy and late-time cosmic expansion using the Weighted Function Regression method

Recent data from multiple supernova catalogs and DESI, when combined with CMB, suggest a non-trivial evolution of dark energy (DE) at the $2.5-4\sigma$ CL. This evidence is typically quantified using the CPL parametrization of the DE equation-of-state parameter which corresponds to a first-order Taylor expansion around $a = 1$. However, this truncation is to some extent arbitrary and may bias our interpretation of the data, potentially leading us to mistake spurious features of the best-fit CPL model for genuine physical properties of DE. In this work, we apply the Weighted Function Regression (WFR) method to eliminate the subjectivity associated with the choice of truncation order. We assign Bayesian weights to the various orders and compute weighted posterior distributions of the quantities of interest. Using this model-agnostic approach, we reconstruct some of the most relevant background quantities, examining the robustness of our results against variations in the CMB and SNIa likelihoods. Furthermore, we extend our analysis by allowing for negative DE. Our results corroborate previous indications of dynamical DE, now confirmed for the first time using the WFR method. The combined analysis of CMB, BAO, and SNIa data favors a DE component that transitions from phantom to quintessence at redshift $z_{\rm cross}\sim 0.4$. The probability of phantom crossing lies between 96.21% and 99.97%, depending on the SNIa data set used, and hence a simple monotonic evolution of the DE density is excluded at the $\sim 2-4\sigma$ CL. Moreover, we find no significant evidence for a negative dark energy density below $z\sim 2.5-3$. Our reconstructions do not address the Hubble tension, yielding values of $H_0$ consistent with the Planck/$\Lambda$CDM range. If SH0ES measurements are not affected by systematic biases, the evidence for dynamical dark energy may need to be reassessed. [abridged]