Pseudodifferential equations on the unit sphere in $\R^n$, $n\ge 3$, are considered. The class of pseudodiffrential operators have long been used as a modern and powerful tool to tackle linear boundary-value problems. These equations arise in geophysics, where the sphere of interest is the earth. Efficient solutions to these equations on the sphere become more demanding when given data are collected by satellites. In this dissertation, firstly we solve these equations by using spherical radial basis functions. The use of these functions results in meshless methods, which have recently become more and more popular. In this dissertation, the collocation and Galerkin methods are used to solve pseudodifferential equations. From the point of view of application, the collocation method is easier to implement, in particular when the given data are scattered. However, it is well-known that the collocation methods in general elicit a complicated error analysis. A salient feature of our work is that error estimates for collocation methods are obtained as a by-product of the analysis for the Galerkin method. This unified error analysis is thanks to an observation that the collocation equation can be viewed as a Galerkin equation, due to the reproducing kernel property of the space in use. Secondly, we solve these equations by using spherical splines with Galerkin methods. Our main result is an optimal convergence rate of the approximation. The key of the analysis is the approximation property of spherical splines as a subset of Sobolev spaces. Since the pseudodifferential operators to be studied can be of any order, it is necessary to obtain an approximation property in Sobolev norms of any real order, negative and positive. Solving pseudodifferential equations by using Galerkin methods with spherical splines results, in general, in ill-conditioned matrix equations. To tackle this ill-conditionedness arising when solving two special pseudodifferential equations, the Laplace--Beltrami and hypersingular integral equations, we solve them by using a preconditioner which is defined by using the additive Schwarz method. Bounds for condition numbers of the preconditioned systems are established.