When the rate variation among sites is described by a gamma distribution, an important problem is how to estimate the shape parameter alpha, which is an index of the degree of among-site rate variation. The parsimony-based methods for estimating alpha are simple but biased, i.e., alpha tends to be overestimated. On the other hand, the likelihood-based methods are asymptotically unbiased but take a huge amount of computational time. In this paper, we have developed a new method to solve this problem: we first estimate the expected number of substitutions at each site, which is corrected for multiple hits, and then estimate the parameter alpha. Our method is computationally as fast as the parsimony method, and the estimation accuracy is much higher than that of parsimony and similar to that of the likelihood method.