Problems and solutions for estimating indel rates and length distributions

Abstract

Insertions and deletions (indels) are fundamental but understudied components of molecular evolution. Here we present an expectation-maximization algorithm built on a pair hidden Markov model that is able to properly handle indels in neutrally evolving DNA sequences. From a data set of orthologous introns, we estimate relative rates and length distributions of indels among primates and rodents. This technique has the advantage of potentially handling large genomic data sets. We find that a zeta power-law model of indel lengths provides a much better fit than the traditional geometric model and that indel processes are conserved between our taxa. The estimated relative rates are about 12-16 indels per 100 substitutions, and the estimated power-law magnitudes are about 1.6-1.7. More significantly, we find that using the traditional geometric/affine model of indel lengths introduces artifacts into evolutionary analysis, casting doubt on studies of the evolution and diversity of indel formation using traditional models and invalidating measures of species divergence that include indel lengths.

FIG. 1.—

SMUVE, a generalized pair HMM architecture. Match state (M) emits a pair of nucleotides to the alignment via the K2P model, whereas indel states (U and V) emit nucleotides to only one sequence based on either a zeta or geometric distribution. Transition probabilities are evenly distributed out of each node unless otherwise noted; τ = 1/(a + 1) and η = e^−2rt.

FIG. 2.—

Comparison of optimal and EM estimation procedures. Each line represents one of four procedures: geo-SMUVE, geometric optimal, zeta-SMUVE, and zeta-optimal. RMSE and biases have been normalized using the true value of the parameters. Note that the indel length distribution parameters of the geometric and zeta models are not directly comparable to one another.

FIG. 3.—

Estimated maximum likelihood SMUVE parameters, with 95% confidence intervals. Each bar is an estimate from a pairwise comparison among human, chimp, mouse, or rat. Distribution parameters are not directly comparable between the two models and are plotted on different axes.

FIG. 4.—

Histograms of the observed number of indels of lengths 1–9 in the mouse–rat alignment space, with 95% confidence intervals. Triangles represent MLEs.