A modified form of Shepard interpolation of ab initio molecular potential energy surfaces is presented. This approach yields significant improvement in accuracy over previous related schemes. Here each Taylor expansion used in the interpolation formula is assigned a confidence volume which controls the relative weight assigned to that expansion. The parameters determining this confidence volume are derived automatically from a simple Bayesian analysis of the interpolation data. As the iterative scheme expands the data set, the confidence volumes are also iteratively refined. The potential energy surfaces for nine reactions are used to illustrate the accuracy obtained.