  We alleviate the circularity problem, whereby gamma-ray bursts are not
perfect distance indicators, by means of a new model-independent technique
based on B\'ezier polynomials. To do so, we use the well consolidate
\textit{Amati} and \textit{Combo} correlations. We consider improved calibrated
catalogs of mock data from differential Hubble rate points. To get our mock
data, we use those machine learning scenarios that well adapt to gamma ray
bursts, discussing in detail how we handle small amounts of data from our
machine learning techniques. In particular, we explore only three machine
learning treatments, i.e. \emph{linear regression}, \emph{neural network} and
\emph{random forest}, emphasizing quantitative statistical motivations behind
these choices. Our calibration strategy consists in taking Hubble's data,
creating the mock compilation using machine learning and calibrating the
aforementioned correlations through B\'ezier polynomials with a standard
chi-square analysis first and then by means of a hierarchical Bayesian
regression procedure. The corresponding catalogs, built up from the two
correlations, have been used to constrain dark energy scenarios. We thus employ
Markov Chain Monte Carlo numerical analyses based on the most recent Pantheon
supernova data, baryonic acoustic oscillations and our gamma ray burst data. We
test the standard $\Lambda$CDM model and the Chevallier-Polarski-Linder
parametrization. We discuss the recent $H_0$ tension in view of our results.
Moreover, we highlight a further severe tension over $\Omega_m$ and we conclude
that a slight evolving dark energy model is possible.
