Schuster's two-stream approximation is first derived from Chandrasekhar's radiative transfer equation and is then generalized to an arbitrary number of streams. The resulting technique for solving the transfer equation which is similar to the discrete ordinate and spherical harmonic methods, is found to be especially useful for modelling atmospheres with complicated phase functions and moderate optical depths. To illustrate the method, a four-stream approximation is evaluated for a Henyey-Greenstein phase function with an asymmetry factor of g = 0.5.