Nearest-neighbor distribution functions characterize the probability of finding a nearest neighbor at some given distance from a reference point in systems of interacting particles and are of fundamental importance in a variety of problems in the physical and biological sciences. We extend the formalism of Torquato, Lu, and Rubinstein [Phys. Rev. A 41, 2059 (1990)] for identical spheres to obtain exact series representation of nearest-neighbor functions (void and particle probability densities) and closely related quantities for systems of interacting D-dimensional spheres with a polydispersivity in size. Polydispersivity constitutes a basic feature of the structure of random systems of particles and leads to a wider choice of possible definitions for nearest-neighbor functions. The most relevant definition for a polydispersed system of particles involves the ``nearest particle surface'' rather than the ``nearest particle center'' and thus we refer to them as ``nearest-surface distribution functions.'' For the special cases of D-dimensional hard and overlapping spheres, we obtain analytical expressions for the nearest-surface functions that are accurate for a wide range of sphere concentrations. Employing these results, we are able to compute the corresponding mean nearest-surface distances for polydispersed hard spheres. Finally, we determine the nearest-surface functions for bidispersed systems from Monte Carlo computer simulations and find that our theoretical results are in very good agreement with the data.