where K y = K + σ 2 y I the covariance matrix with some noise and y = f + ϵ the training dataset with also noise defined as ϵ ∼ N (0, σ 2 y I) which follows a Gaussian distribution, and added independently to each observation. Unfortunately, when applying GPR models on a dataset of size n, the exact inference has complexity O(n 3 ) with storage demands of O(n 2 ).