Context. We present a probabilistic cross-correlation approach to estimate time delays in the context of reverberation mapping (RM) of Active Galactic Nuclei (AGN).

Aims. We reformulate the traditional interpolated cross-correlation method as a statistically principled model that delivers a posterior distribution for the delay.

Methods. The method employs Gaussian processes as a model for observed AGN light curves. We describe the mathematical formalism and demonstrate the new approach using both simulated light curves and available RM observations.

Results. The proposed method delivers a posterior distribution for the delay that accounts for observational noise and the non-uniform sampling of the light curves. This feature allow us to fully quantify its uncertainty and propagate it to subsequent calculations of dependant physical quantities, e.g., black hole masses. It delivers out-of-sample predictions, which enables us to subject it to model selection and it can calculate the joint posterior delay for more than two light curves.

Conclusions. Because of the numerous advantages of our reformulation and the simplicity of its application, we anticipate that our method will find favour not only in the specialised community of RM, but in all fields where cross-correlation analysis is performed. We provide the algorithms and examples of their application as part of our Julia GPCC package.