We present an entirely new class of high-order numerical algorithms for computational fluid dynamics. The new method is based on the Gaussian Processes (GP) modeling that generalizes the Gaussian probability distribution. Our approach is to adapt the idea of the GP prediction technique which utilizes the covariance kernel functions, and use it to reconstruct a high-order approximations for computational simulations. We propose the GP high-order method as a new numerical high-order formulation, alternative to the conventional polynomial-based approaches. We will show that the GP method is shown to be much faster in both convergence and performance rates than the popular choices of polynomial-based high-order methods such as PPM, WENO-5, and WENO-Z.