The current paper introduced two approximation operators of large scattered datasets for spherical interpolation. The suggested solution method is an extension of Shepard's well-known method of spherical interpolating, which uses the inverted distances of scattered points as weight functions. With regard to this, the first proposed operator is a linear combination of basis functions with coefficients that are the values of the function. As for the second operator, we consider a spherical triangulation of the scattered points and substitute function values with a local interpolant, which locally interpolates the given data at the vertices of each triangle. Moreover, numerical tests have been carried out to show the interpolation performance, where several numerical results reveal the signified approximation accuracy of the proposed operators.