## Abstract

The allometric-constraint hypothesis states that evolutionary divergence of morphological traits is restricted by integrated growth regulation. In this study, we test this hypothesis on a time-calibrated and well-documented palaeontological sequence of dental measurements on the Pleistocene arvicoline rodent species *Mimomys savini* from the Iberian Peninsula. Based on 507 specimens representing nine populations regularly spaced over 600 000 years, we compare static (within-population) and evolutionary (among-population) allometric slopes between the width and the length of the first lower molar. We find that the static allometric slope remains evolutionary stable and predicts the evolutionary allometry quite well. These results support the hypothesis that the macroevolutionary divergence of molar traits is constrained by static allometric relationships.

## 1. Introduction

Huxley's [1,2] model of allometric growth provides a unifying theoretical framework to link within-species variational properties to morphological divergence. Based on the assumption that the growth rates of morphological traits are proportional to the growth of body size, this model predicts that the trait size, *Y*, will covary in the population with body size, *X*, according to a power relationship, *Y* = *aX ^{b}*, that becomes linear on a log scale: log[

*Y*] = log[

*a*] +

*b*log[

*X*] [2]. The value of the allometric slope

*b*expresses the growth rate of the trait relative to the growth rate of body size [3–6].

The allometric scaling model occupies a key place in studies of developmental constraints on morphological evolution (e.g. [3,4]). Originally based on a developmental model, morphological allometry can be defined not only at the ontogenetic level, but also at the static and evolutionary levels, that is, within and among evolutionary units (i.e. populations, species or higher taxa), respectively. The allometric-constraint hypothesis postulates that the static allometric slope (*b*) remains stable at macroevolutionary timescales and is able to restrict the evolutionary divergence in morphospace along its specific trajectory (for reviews, see [5,6]). This narrow-sense static allometric slope may therefore be a specific case of a ‘genetic line of least resistance’ (*sensu* Schluter [7]) explicitly encompassing developmental or physiological constraints.

Nowadays, this hypothesis is often expressed in a quantitative-genetics framework, and it has been predicted that bivariate drift or selection on body size alone can generate an evolutionary allometric slope collinear with the genetic static allometric slope [8,9]. Recently, comparative studies have provided substantial arguments for the stability and the constraining role of the narrow-sense static allometric slope [6,10,11]. Few studies, however, have looked at the evolution of ontogenetic or static allometric slopes in the fossil record.

Although evolutionary sequences in the fossil record often display a stationary pattern [4,12,13], there are also many examples of lineages exhibiting apparent ‘trends’ of size change. These have been attributed to constant directional selection (e.g. [14]), undirected diversification from a lower size bound [15] or to release of constraints associated with the colonization of new ecological niches (e.g. [16]). Examples of bivariate evolutionary divergence of populations or species distributed along straight lines (on arithmetic or logarithmic scales) have been documented from early studies of fossil sequences (e.g. [17–20]) with the allometric constraint as a background hypothesis. However, the difficulty of estimating the static allometric slope from limited fossil samples, as well as the lack of adequacy of the approaches generally used [5, 6], have hampered evaluations of the allometric-constraint hypothesis from palaeontological sequences. Filling this gap requires accounting for time dependency of fossil samples and a dataset large enough to properly estimate both static and evolutionary allometries in order to compare them.

The molars of arvicoline rodents are well suited for the study of evolutionary trends in the fossil record. They are usually well-preserved fossils found in abundance in Neogene fossiliferous layers and they are complex and taxonomically diagnostic organs for which a large number of traits can be measured and compared (e.g. [21–24]). In the following, we test the allometric-constraint hypothesis by estimating how the static allometry predicts the evolutionary allometry between molar traits in a recently described 600 000 year-long fossil sequence for a Pleistocene arvicoline rodent species in the Iberian Peninsula [25,26].

## 2. Material and methods

### (a) Samples

We reanalysed data from Lozano-Fernández *et al.* [25,26] describing 600 000 years of microevolutionary change of the first lower molar (M_{1}) dimensions of the rodent species *Mimomys savini*. This extint arvicoline is considered as the first representative of the water-vole lineage (genus *Arvicola*). *Mimomys savini* occurred in Europe during the Early Pleistocene and beginning of the Middle Pleistocene, between 1.8 and 0.6 Ma [27–31]. In the Iberian Peninsula, the oldest record of this species is located in the deposits of the Guadix-Baza basin [32–34].

*Mimomys savini* exhibits evolutionary trends, including an increase in size, in both central Europe and the Iberian Peninsula [25,31,35–38]. However, individuals in the Iberian populations are systematically larger, as supported by the comparison of means from Lozano-Fernández *et al.* [25,26] with means from Maul *et al.* [31,36,37], suggesting separate evolutionary histories between the Iberian and Central-European lineages as expected from their geographical isolation [26]. Among Iberian populations, the tendency towards size increase was documented for both the length and the width of the first lower molar [25], suggesting a length–width evolutionary allometric relationship that motivated this study.

The fossil samples come from the Lower Pleistocene layers of different localities in Spain. As their ages coincide with early human dispersal events, the Iberian *M. savini* samples have indirectly benefited from absolute dating efforts providing an exceptionally long and high-quality chronology for a terrestrial Pleistocene lineage. The sequence consists of nine samples associated with absolute ages. Their chronology (table 1) was taken from Lozano-Fernández *et al.* [25,26] and updated according to recent absolute dating [39,40]. Molar length and width (figure 1) have been measured (in millimetres) for a total of 507 individuals as described by Lozano-Fernández *et al.* [25].

The number of sampled individuals per population ranges from nine to 137 (median: 51). Because five out of the nine samples have a size exceeding 50 individuals, this dataset is well suited to estimate static allometric slopes, and to assess their variation through time and their effect on evolutionary allometry. Some of the within-population variation in molar morphology may be due to variation in dental wear. For example, the complexity of the M_{1} occlusal outline has been shown to correlate with wear in arvicoline rodents [41]. In this study, however, extreme stages of wear were removed from the sample, and the parallel crown walls of the *M. savini* hypselodont M_{1} make wear unlikely to have a large influence on the length–width pattern of static allometry.

Here, we analysed the allometric relationship between the two traits, regressing the natural log width of the first molar on its natural log length.

### (b) Statistics

We first estimated variation in static allometric slopes among populations with a mixed-effects regression model, fitted in the R package *lme4* [42], with log of molar width as response variable, log of molar length as predictor variable and population as a random factor. Two models fitted with restricted maximum likelihood (REML) were compared using the Akaike's information criterion (AIC), one assuming a static allometric slope common for all populations and the other assuming a variable allometric slope across populations. This comparison allowed estimating the among-population differences in static allometric slopes. Because the common-slope model showed the best fit (see Results), we compared the evolutionary allometric slope to the common static slope. To estimate the evolutionary and the common static allometric slope, while accounting for phylogenetic relatedness among populations and measurement error in the population means, we fitted a bivariate Brownian-motion model of evolution using the *MCMCglmm* R package [43]. This model estimates the variance matrix between the two traits, log(width) and log(length), at both the among- and within-population levels. The evolutionary allometric slope (*b*_{e}) was estimated as the ratio of the traits' covariance to the variance in log molar length, based on the phylogenetically corrected among-population variance matrix. The common static allometric slope (*b*_{s}) was equivalently estimated from the within-population variance matrix. From this model, we also estimated the contrast between the slopes at the two levels (i.e. *b*_{e} – *b*_{s}) to test the allometric-constraint hypothesis. The model ran for 250 000 iterations including a burn-in of 50 000 iterations and a thinning interval of 200 iterations to minimize the level of autocorrelation in the Markov chain Monte Carlo (MCMC) resampling. For the variance matrices, we used scaled *F*-distributed priors, where the variances divided by 10 are *F*_{1,1} distributed. For the residuals, we used inverse-Wishart distributed priors with parameters **V** and *n* [43]. The matrix parameter **V** was set to a crude guess of the within-population variance matrix, while *n* was set to 1.002. Confidence intervals for the model parameters were obtained by MCMC resampling.

The expected among-population covariance structure was computed from a phylogeny of populations (handled with the *ape* R package [44]) separated along the stem by internal branches whose lengths correspond to the time elapsed between them (table 1 and figure 2*a*). To account for the possibility that the observed populations may not stand in a perfect ancestor–descendant relation to each other, or that an unquantified factor (e.g. geography) has influenced the ancestor–descendent divergence pattern, we added a non-zero branch length from the stem to each population. The length of this branch was treated as a parameter to be estimated. Hence, we assumed that the divergence from the stem was the same for all populations (except the terminal population) to keep the spacing of populations constant in time. To estimate the length of the branch, we assessed 50 values between 1000 and 100 000 years. Each model was run five times for each value. The best fit according to the deviance information criterion [45] (DIC, averaged over the five replicates) was a terminal branch length of 18 000 years (figure 2*b*). This was retained for the final analysis.

We evaluated whether the pattern of size increase in the molar measurements [25] can be attributed to a directional evolutionary process using Hunt's [46] maximum-likelihood framework for comparing models fitted to fossil time series. We fitted two models: (i) an undirectional Brownian motion and (ii) a directional Brownian motion including a trend parameter. The models were fitted to both log(length) and log(width) independently and compared for each trait with small-sample-corrected AIC values (AICc) in the *paleoTS* R package [46]. For these analyses, we used a model parametrization based on the joint distribution of mean trait values, because it performs better for detecting directional trends [47].

## 3. Results

The least-squares estimates of the static allometric slopes varied across *M. savini* populations (table 1 and figure 3*a*). For the seven populations with a decent sample size (*N* > 20), the static slope ranged between 0.50 and 1.02, with a median of 0.58 (table 1; electronic supplementary material, figure S1). This variation may be due to estimation error, however, because the common-slope model fitted the data better that the variable-slope model (ΔAIC = 2.82).

The common static allometric slope, *b*_{s} = 0.57 (95% CI: 0.47–0.67), was close to the phylogenetically controlled evolutionary allometric slope *b*_{e} = 0.49 (95% CI: −0.11–1.16; table 2). This similarity can be expressed as a small and statistically non-significant contrast between the two slopes (*b*_{e} − *b*_{s} = −0.08, 95% CI: −0.83–0.49; table 2). Figure 3*b* illustrates this pattern, showing that the directions of stepwise inter-population divergence through time in the bivariate plane roughly follow the common static allometric line, leading to an evolutionary allometry nearly collinear with the static allometry.

Between the first and the last point of the *M. savini* fossil sequence, the length and width of the first molar increased by 12% and 7%, respectively (table 1 and figure 3). For both molar length and width, the parameters of directional trends (*μ*_{trait} ± s.e.) are small (*μ*_{log(length)} = 0.15 ± 0.08; *μ*_{log(width)} = 0.35 ± 0.32). Accordingly, the fit of a non-directional Brownian-motion model received a better support with lower AICc values for both molar length (ΔAICc = 2.00) and width (ΔAICc = 3.71). We also tested two other evolutionary models that showed worse fit. Relative to the Brownian motion, an Ornstein–Uhlenbeck model with a single optimum was worse with ΔAICc of 9.02 and 10.29 for M_{1} log(length) and log(width), respectively. Similarly, a white-noise model (‘stasis’) showed a worse fit with ΔAICc of 9.95 and 13.03, respectively.

## 4. Discussion

Comparing narrow-sense static allometric slopes within fossil populations with the evolutionary slope generated during population divergence over time is an intuitive way to test the allometric-constraint hypothesis. In the case of *M. savini*, we found that the static allometric slope between the length and the width of the first lower molar remains stable through time and predicts the evolutionary change observed over 600 000 years. Thus, the size increases previously described for both molar traits by Lozano-Fernández *et al.* [25,26] could be the result of morphological divergence constrained along a static allometric line. Despite its large standard error, the evolutionary allometry is very close to the static allometry, indicating that the allometric-constraint hypothesis is not falsified in this case. This is in accordance with the conclusions from the recent meta-analysis by Voje *et al.* [6] supporting a constraining role of static allometry on trait divergence on timescales less than a million years. Over larger time and taxonomic scales, larger changes in allometric coefficients are commonly observed [6], as claimed for rodent skull traits [48,49].

In figure 4, we reanalysed in a narrow-sense allometric context the data from Lich [24], describing a *ca* 119 000 year-long fossil sequence for the New World rodent species *Cosomys primus*, ‘a morphologically convergent form of the Eurasian *Mimomys*’ [38]. We found a strong M_{1} length–width evolutionary allometry (ordinary least-squares slope: *b*_{e} = 0.94 ± 0.15) close to the one estimated in *M. savini* (*b*_{e} = 0.81 ± 0.18, figure 3*a*). That an analogous pattern of relative trait divergence is observed independently in two distinct fossil lineages adds additional support to a role for allometric constraints in the evolution of molar proportions in rodents.

Evaluating previous work on allometric divergence in rodent molars is challenging because of the diversity of approaches, focal traits and definitions of allometry used by different researchers. Many discussions based on fossil time series have involved contrasting patterns between ‘size’ and ‘shape’, where size has been argued to be more evolutionary labile than shape (e.g. [12,14,50–52]). Such comparison can lead to meaningless conclusions if size and shape are expressed on different scales. In an elementary bivariate case as in the *M. savini* molar, shape is often assessed as a dimensionless ratio (e.g. width/length). Variation in such a ratio depends on the allometric relationship linking the two traits, since *Y*/*X* = *aX ^{b}*

^{−1}[53], and the possibility for shape divergence among populations then depends on how much the evolutionary allometric slope deviates from isometry (i.e. from

*b*= 1). In spite of a long history of critique [54], analyses focusing on ratios and their comparison to size generally ignore their intrinsically underlying allometry. To become meaningful, such comparisons would require standardization to a common scale, but how to do this is not obvious. We suggest it is better to work within the (narrow-sense) allometric framework to understand how the whole size–shape relationship is evolving.

Results from some previous studies indicate that allometry is a strong component of shape variation in rodent molars. For instance, in a fossil series of the arvicoline *Terricola savii*, the length-to-width ratio increased over time [23] and also varies among contemporary populations [55] indicating a length–width evolutionary allometry deviating from isometry (i.e. *b* ≠ 1). In *Myodes glareolus* [56], another arvicoline, and within the *Arvicola* genus [57], a strong effect of size on M_{1} shape suggests an analogous pattern. Other rodent taxa exhibiting a molar pattern convergent with arvicoline rodents could be informative. For instance, the spectacular elongation of the lower first molar found in the *Mikrotia* murine genus endemic to the Gargano palaeo-archipelago [58] could similarly be explained by a hypo-allometric constraint (i.e. *b* < 1) expressed during the size increase documented in some of these lineages [22,58] or result from selection imposed by a specific food regime. Evaluating allometric constraints in this taxon could shed light on how such an exaggerated dental morphology evolves.

With the caveat that the *M. savini* time series has a small number of sampling intervals and the power to detect evolutionary patterns is consequently low, we did not detect a trend or any more complex patterns in molar size evolution. Our tentative interpretation is thus that the molars diversify according to a Brownian motion along a fixed allometric relationship. We note, however, that the static allometries were far from tight. For the better sampled populations, the *R*^{2} ranges from 15 to 34% with a *R*^{2} of 20% for the estimated common slope. This does not fit well with the allometric-constraint hypothesis [6]. It is hard to see how a trait can be substantially constrained by another trait that only explains 20% of its variance. After all, 80% of the variance should be freely available for selection to change the traits away from the allometric relation. Part of the low *R*^{2} must be due to measurement error, which is unquantified in these data, but this can surely not explain the big pattern. This leaves us with the conundrum that we seem to find evidence for allometric constraints resulting from a weak allometric relationship. A possible solution to this conundrum may be to view the allometric relations as resulting from selective rather than genetic constraints. If both the evolutionary and the static allometries in *M. savini* have been shaped by similar selective pressures, then they would be similar even without a tight fit of the traits to the static allometry. This ‘adaptationist’ hypothesis posits that the populations are aligned along a selective ridge (*sensu* Armbruster & Schwaegerle [59]) that also imposes stabilizing selection on their static allometric slopes. This scenario matches Arnold *et al.*'s [60,61] concept of a ‘selective line of least resistance’. Thus, the static allometry could not only be a direct predictor of the most likely response to selection but also a local reflection of the dominant features of the adaptive landscape. In both cases, the consequence is that estimates of the static allometric slopes contain information to predict evolutionary divergence.

To summarize, we found that the relative evolutionary increase in the size of two molar traits in *M. savini* is consistently predicted by a static allometric relationship remaining stable over more than half a million year. We hope this result will encourage further investigations on other fossil time series for which allometric slopes can reliably be estimated at several levels.

## Data accessibility

Data available from the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.mg0v4.

## Funding statement

This work was supported by grant 196434/V40 from the Norwegian Research Council to C.P. at NTNU, a pre-doctoral grant from the Fundación Atapuerca assigned to the IPHES (Institut Català de Paleoecologia Humana i Evolució Social) and projects CGL2012-38358, CGL2012-38434-C03-01, BP2012-CGL38358 and CGL 2009-12703-C03-03 of the Spanish Ministry of Economy and Competitiveness, SGR 2009-324 of the Generalitat de Catalunya and the University of Zaragoza, Gobierno de Aragon, European Social Fund.

## Acknowledgements

This Spanish–Norwegian collaboration was initiated within the scientific team of the Bois-de-Riquet archaeological site (dir. L. Bourguignon and J.-Y. Crochet). We thank Gene Hunt and an anonymous referee for constructive comments on the manuscript. The sites’ excavation team has helped with the extraction, sieving and washing of sediments.

## Footnotes

One contribution of 14 to a Theme Issue ‘Phenotypic integration and modularity in plants and animals’.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.