We discuss a Bayesian methodology for the solution of the inverse problem
underlying the determination of parton distribution functions (PDFs). In our
approach, Gaussian Processes (GPs) are used to model the PDF prior, while Bayes
theorem is used in order to determine the posterior distribution of the PDFs
given a set of data. We discuss the general formalism, the Bayesian inference
at the level of both parameters and hyperparameters, and the simplifications
which occur when the observable entering the analysis is linear in the PDF. We
benchmark the new methodology in two simple examples for the determination of a
single PDF flavor from a set of Deep Inelastic Scattering (DIS) data and from a
set of equal-time correlators computed using lattice QCD. We discuss our
results, showing how the proposed methodology allows for a well-defined
statistical interpretation of the different sources of errors entering the PDF
uncertainty, and how results can be validated a posteriori