For the interpolation of graph signals with generalized shifts of a graph basis function (GBF), we introduce the concept of positive definite functions on graphs. This concept merges kernel-based interpolation with spectral theory on graphs and can be regarded as a graph analog of radial basis function interpolation in euclidean spaces or spherical basis functions. We provide several descriptions of positive definite functions on graphs, the most relevant one is a Bochner-type characterization in terms of positive Fourier coefficients. These descriptions allow us to design GBFs and to study GBF interpolation in more detail: we are able to characterize the native spaces of the interpolants, we provide explicit estimates for the interpolation error and obtain bounds for the numerical stability. As a final application, we show how GBF interpolation can be used to get quadrature formulas on graphs.