We present a finite dimensional Gaussian process (GP) regression for no-arbitrage interpolation and modeling the local volatility surface. In this setup, the MAP estimate serves the purpose of identifying the locations of the most likely arbitrages in the data and quantifying them. Furthermore, Hamiltonian Monte Carlo can be used to efficiently sample from the posterior price surface and provide UQ of the local volatility surface. We demonstrate the performance of this approach relative to SSVI and a NN approach on SPX options. This is joint work with Stéphane Crépey and Areski Cousin.