Markov chain Monte Carlo techniques applied to parton distribution functions determination: Proof of concept

We present a new procedure to determine parton distribution functions (PDFs), based on Markov chain Monte Carlo (MCMC) methods. The aim of this paper is to show that we can replace the standard ${\chi }^2$ minimization by procedures grounded on statistical methods, and on Bayesian inference in particular, thus offering additional insight into the rich field of PDFs determination. After a basic introduction to these techniques, we introduce the algorithm we have chosen to implement—namely Hybrid (or Hamiltonian) Monte Carlo. This algorithm, initially developed for Lattice QCD, turns out to be very interesting when applied to PDFs determination by global analyses; we show that it allows us to circumvent the difficulties due to the high dimensionality of the problem, in particular concerning the acceptance. A first feasibility study is performed and presented, which indicates that Markov chain Monte Carlo can successfully be applied to the extraction of PDFs and of their uncertainties.

Figure 1

Values of the parameter $B_g$ as a function of the Monte Carlo time for three independent Markov chains. The green solid line represents a chain starting from the value given by MINUIT minimization, whereas the two other chains start from values much higher (blue dotted) or much lower (red dashed). The initial corresponding ${\chi }^2/\mathrm{d}.\mathrm{o}.\mathrm{f}.$ values are, from chain 1 to chain 3, respectively, ${\chi }^2/\mathrm{d}.\mathrm{o}.\mathrm{f}.=67.44$, ${\chi }^2/\mathrm{d}.\mathrm{o}.\mathrm{f}.=0.87$ and ${\chi }^2/\mathrm{d}.\mathrm{o}.\mathrm{f}.=81.58$. We identify clearly on this plot the thermalization region, which is limited to the first $\sim 100–210$ iterations.

Figure 2

Values of the parameter $B_g$ as a function of the Monte Carlo time for three independent Markov chains (lhs). The starting points of the chains and the color code are the same as in Fig. 1. The chains are clearly converging towards the same stationary distribution, what is confirmed by plotting the parameter distribution for each chain (rhs), after removing thermalization points and taking into account autocorrelation.

Figure 3

Probability distribution functions of the PDF parameters (diagonal) and 2D correlation plots between parameters (off-diagonal).

Figure 4

Marginal probability distribution of parameters $B_g$ (lhs) and $C_{\overline{D}}$ (rhs). They do not follow a Gaussian law, as can be seen from the gaussian fit (solid red line).

Figure 5

${\chi }^2$ distribution for a 10D Monte Carlo chain. The solid red line is an adjustment of these data with a ${\chi }^2$ distribution law with 10 degrees of freedom. The dashed (dashdot) vertical line indicates the 68% (95%) confidence limit.

Figure 6

Gluon PDF probability distribution function for $x\approx {10}^{-4}$ (lhs) and $x\approx 0.83$ (rhs) at fixed $Q^2=10{\mathrm{GeV}}^2$. The 68% confidence interval is obtained from this distribution, considering the region of the distribution containing 68% of the data remaining on each side of the best fit value.

Figure 7

The parton distribution functions obtained using MCMC (right) compared to HERAPDF1.0 (ZMVFN scheme) from xFitter output (left) for ${\mathrm{xu}}_{\mathrm{val}}$ and xg, at ${\mathrm{Q}}^2=10{\mathrm{GeV}}^2$. The bands show the 68% confidence interval around the central value (in solid red line) for the MCMC PDFs, and the standard $\mathrm{\Delta }{\chi }^2=1$ deviation for HERAPDF.

Figure 8

Ratio of MCMC PDFs and HERAPDF1.0 (ZMVFN scheme) central values for ${\mathrm{xu}}_{\mathrm{val}}$ and xg at ${\mathrm{Q}}^2=10{\mathrm{GeV}}^2$.

Figure 9

Comparison of the PDF uncertainties, normalized by the best fit value, as determined by the Hessian and MCMC methods at NLO for the valence distribution $x{u}_{\mathrm{val}}$ and the gluon distribution $xg$, at a scale ${\mathrm{Q}}^2=10{\mathrm{GeV}}^2$.

Figure 10

The MCMC parton distribution functions $\mathrm{x}\overline{u}$, $\mathrm{x}\overline{d}$, $\mathrm{x}\overline{s}$ and $\mathrm{x}\overline{c}$ at ${\mathrm{Q}}^2=10{\mathrm{GeV}}^2$. The bands show the 68% confidence interval around the best fit value (in solid red line).

Figure 11

PDFs of valence quarks ($x{u}_{\mathrm{val}}$, $x{d}_{\mathrm{val}}$), sea quarks ($x\mathrm{\Sigma }=x\overline{u}+x\overline{d}+x\overline{s}+x\overline{c}$, with $x\overline{c}=0$ for $Q^2<m_c^2$), and gluon ($xg$) using MCMC at the scale ${\mathrm{Q}}^2=1.9{\mathrm{GeV}}^2$ (lhs) and ${\mathrm{Q}}^2=10{\mathrm{GeV}}^2$ (rhs) obtained with the ZMVFN scheme. The gluon and the sea distributions are scaled down by a factor of 20. The experimental uncertainties (68% confidence limit as defined in the text from the probability density of PDFs) are represented by the green-shaded region.