doi: 10.1111/jeb.12523
Evolution of wing shape in hornets: why is the wing venation efficient for species identification? A. PERRARD*†, M. BAYLAC†, J. M. CARPENTER* & C. VILLEMANT† *Division of Invertebrate Zoology, American Museum of Natural History, New York, NY, USA †UMR 7205 ISYEB, Museum National d’Histoire Naturelle, Paris, France
Keywords:
Abstract
Gaussian mixture models; genus Vespa; insect wing; phylogenetic signal; species recognition; tangent shape space.
Wing venation has long been used for insect identification. Lately, the characterization of venation shape using geometric morphometrics has further improved the potential of using the wing for insect identification. However, external factors inducing variation in wing shape could obscure specific differences, preventing accurate discrimination of species in heterogeneous samples. Here, we show that interspecific difference is the main source of wing shape variation within social wasps. We found that a naive clustering of wing shape data from taxonomically and geographically heterogeneous samples of workers returned groups congruent with species. We also confirmed that individuals can be reliably attributed to their genus, species and populations on the basis of their wing shape. Our results suggested that the shape variation reflects the evolutionary history with a potential influence of other factors such as body shape, climate and mimicry selective pressures. However, the high dimensionality of wing shape variation may have prevented absolute convergences between the different species. Wing venation shape is thus a taxonomically relevant marker combining the accuracy of quantitative characters with the specificity required for identification criteria. This marker may also highlight adaptive processes that could help understand the wing’s influence on insect flight.
Introduction Species identification and classification is a first and crucial step in many biological studies. It allows researchers to attribute information to recognizable entities (Mayr, 1969; Gaston & O’Neill, 2004). However, species recognition is often challenged by taxonomic diversity and the subtlety of interspecific variation (Meier et al., 2006; Bickford et al., 2007). With hundreds of thousands of distinct species, insects are one of the most diverse taxa of the living world (Grimaldi & Engel, 2005). Wing venation patterns are classically used for their identification at the order or family levels (Comstock & Needham, 1899). More recently, the quantification of the wing shape through geometric morphometrics became widely Correspondence: Adrien Perrard, Division of Invertebrate Zoology, American Museum of Natural History, 79th Street at Central Park West, New York, NY 10024, USA. Tel.: +1 212 496 3446; fax: +1 212 769 5277; e-mail: [email protected]
used to study subtle variation among insect genera, species and even populations, sex and castes within the same species (e.g. Baylac & Daufresne, 1996; Baylac et al., 2003; Villemant et al., 2007; Bai et al., 2011; De Meulemeester et al., 2012; Perrard et al., 2012; Outomuro et al., 2013). These previous studies showed that wing shape differs between most taxonomic groups and can be used for insect identification (Matias et al., 2001; Steinhage et al., 2001; Francoy et al., 2008, 2009; Henry et al., 2010). However, they did not address how the wing shape varies along a lineage: Is this marker pertinent for a generalized identification, independent from the taxonomic scale of the sampling? Is the morphological divergence between taxa the main source of variation in the venation shape? Are there other factors influencing its evolution? Geometric morphometrics aim at quantifying geometric shapes, a spatial description of objects excluding size, orientation and position (Bookstein, 1991; Adams et al., 2004). Insect wing shapes are usually described by
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
1
2
A. PERRARD ET AL.
landmarks placed at the vein intersections (Baylac & Daufresne, 1996; Klingenberg et al., 2002). Superimposition methods enable comparison of different shapes by extracting a standardized shape description from the raw landmark coordinates. The space of all possible shapes described by the same number of landmarks is a nonlinear multidimensional manifold called the ‘shape space’ (Kendall, 1977). Because most multivariate statistical analyses are designed for linear spaces, the aligned coordinates resulting from the superimposition are usually projected onto a linear tangent space. In biology, the difference between the shape space and the tangent space is considered negligible: biological variation in shape being very limited compared to the geometric possible shapes, the subspace of biological shapes remains close to the tangent space (Marcus et al., 2000; Zelditch et al., 2004). An ensemble of shapes can thus be represented as a cloud of points in the multivariate tangent space. To be a generalized marker for taxon identification, the variation in the wing venation should enable identification at different levels of the classification. It should be highly structured by evolution and should have undergone low convergence between taxa: the wing venation shape should thus be more influenced by the genetic drift between taxa than by other factors of variation. This implies several predictions as to the distribution of wing shapes within the tangent space: first, wing shapes should be distinct between taxa at different taxonomic levels; second, for a given taxonomic level, the intragroup variation should be lower than the intergroup variation. Finally, as a corollary of the two previous statements, the shape difference between distantly related groups (e.g. different genera)
is expected to be higher than between more closely related taxa. The average shape difference between taxa should thus reflect their classification level and, indirectly, their degree of relatedness. To test these hypotheses, we analysed the evolution of the wing shape of social wasps of the subfamily Vespinae at different classification levels, from genera to populations. Wing shape is influenced by external factors affecting the development of organisms (Debat et al., 2006), and this variation could blur the difference between groups, especially in case of phenotypic convergence. To address the full variability of shapes within each group and to test the difference between groups of the same taxonomic level, we used wild specimens collected from different localities for each species or populations. We first explored the structure of the variability of wing shapes in the tangent space using Gaussian-based clustering. Then, we estimated the difference between the taxonomic groups using linear models and cross-validation attribution tests. Finally, we tested the phylogenetic signal of the wing shape to establish whether the variation in wing shape between the taxonomic groups reflected their relatedness.
Materials and methods Social wasp specimens were gathered from natural history collections of eight public institutions and three private collections. We sampled 1025 specimens among the four genera of the subfamily Vespinae with a special emphasis on hornets by including the 22 species of the genus Vespa and nine populations of the yellow-legged hornet Vespa velutina Lepeletier, 1836 (Table 1).
Species
N
Species [Population]
N
Dolichovespula media (Retzius, 1783) Provespa anomala (de Saussure, 1854) Provespa barthelemyi (du Buysson, 1905) Provespa nocturna v. d. Vecht, 1935 Vespa affinis (Linnaeus, 1764) Vespa analis Fabricius, 1775 Vespa basalis Smith, 1852 Vespa bellicosa de Saussure, 1854 Vespa bicolor Fabricius, 1787 Vespa binghami du Buysson, 1905 Vespa crabro Linnaeus, 1758 Vespa ducalis Smith, 1852 , 1884 Vespa dybowskii Andre Vespa fervida Smith, 1858 Vespa fumida v. d. Vecht, 1959 Vespa luctuosa de Saussure, 1854 Vespa mandarinia Smith, 1852 Vespa mocsaryana du Buysson, 1905
10 51 9 11 33 34 27 32 25 24 49 33 21 23 13 19 30 32
rez, 1910 Vespa multimaculata Pe Vespa orientalis Linnaeus, 1771 Vespa philippinensis de Saussure, 1854 Vespa simillima Smith, 1868 Vespa soror du Buysson, 1905 Vespa tropica (Linnaeus, 1758) Vespa velutina Lepeletier, 1836 [Java] V. velutina [Lombok] V. velutina [Flores] V. velutina [Nepal] V. velutina [Sulawesi] V. velutina [Thailand] V. velutina [Vietnam] V. velutina [Yunnan] V. velutina [Zhejiang] Vespa vivax Smith, 1870 Vespula germanica (Fabricius, 1793) Vespula vulgaris (Linnaeus, 1758)
34 38 10 30 34 31 54 40 25 68 10 33 10 35 63 6 15 13
Table 1 Sampling of the species and populations of Vespinae.
Column ‘N’ indicates the sample size.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Right forewings of specimens were measured using the protocol for pinned specimens described in Perrard et al. (2012). Only workers were included to avoid caste dimorphism. Forewing shape was assessed using 19 landmarks recorded with TPSDig2 software (Rohlf, 2010; Fig. 1). As the tangent space is defined by the superimposition, landmarks of all the specimens were superimposed only once using a generalized Procrustes analysis and thus projected in the same tangent space (Dryden & Mardia, 1998). In this study, wing shapes were described as the scores of the 34 principal components (PC) that have eigenvalues > 0 of a principal component analysis (PCA) of these projections. They are referred to as ‘specimens’ in the text. Variability of wing shape in the tangent space We first explored the variability of the wing shape in the tangent space at three classification levels, from population to genera. Because of the sampling ranging from V. velutina populations to Vespinae genera, the species V. velutina and the genus Vespa were overrepresented in the main PCA compared with the other species and genera. To avoid the bias towards these groups when analysing the variation at the different levels of classification, we performed three additional PCA. The first included the 338 specimens of V. velutina from the nine populations sampled (PCAPop). The second included 50 specimens of V. velutina randomly selected among the 338 specimens and every specimen of the 21 other species of Vespa (PCASp). The third PCA included 50 specimens of Vespa randomly selected from the previous sample and representatives of the three other Vespinae genera: Dolichovespula, Provespa and Vespula (PCAGen). We compared the variances of these subsamples and of the total sample using the respective sums of the 34 PC scores’ variances normalized by the ratio between each subsample size and the total sample size.
3
Gaussian mixture models We looked for groups of similar shapes in the whole sample and in the different subsamples using a Bayesian clustering method called Gaussian mixtures models (GMM; Fraley & Raftery, 1998). This method requires no a priori information on groups or number of groups; it clusters individuals based on the likelihood of groups of Gaussian distributions in the tangent space. If the divergence between taxonomic groups was the main cause of shape variation, Gaussian clusters should be congruent with these groups. Differences between taxa The difference in wing shapes between groups was estimated by a reassignment method based on canonical variate analyses (CVA) with leave-one-out cross-validation. Under the hypothesis of a complete separation of the groups within the tangent space, a perfect attribution would be expected with cross-validated reassignment. Because several groups had small sample sizes compared with the number of variables, the validity of CVA classifications was investigated (Mitteroecker & Bookstein, 2011). Reassignments were performed with different number of variables and excluding progressively the least sampled groups. The results were compared to reassignment obtained from the k-nearest neighbours nonparametric method with k = 1, 3 and 5 (Ripley, 1996). The amount of variation at the different classification levels was determined on the whole sample by the sum of variances (∑2) between the mean shapes of (i) the four Vespinae genera; (ii) the 22 Vespa species; and (iii) the nine populations of V. velutina. Under the hypothesis of a neutral evolution of the shape driven by genetic drift, the sum of squares of mean shapes of higher classification levels should be higher than those corresponding to lower levels or intragroup variation. Phylogenetic signal
Fig. 1 Location of the 2D landmarks measured on the right forewing of specimens.
The congruence between the wing shape variation and the evolutionary history was tested using the phylogenetic signal test from Klingenberg & Gidaszewski (2010). The phylogenetic tree was mapped within the tangent space by computing the ancestral states minimizing the squared parsimony (Huey & Bennett, 1987). The squared length of the tree within the tangent space (D2) was then used to test the phylogenetic signal with a significance level of 0.05. Under a Brownian motion model of evolutionary change of the wing, the phylogeny should be significantly shorter than random trees in the tangent space. The null distribution of the squared tree lengths was computed by permutation of the taxa across the phylogeny.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
c
d
0.02
b
(c)
d 17
14
7
13
27 23 8 24 9 20 21 16 26 10 25 15 22
11
0.02 PC2 (12.6 %) −0.04 −0.02 0.00
3 12 c 6
18
19
−0.02 −0.01 0.00 0.01 0.02 PC1 (30.35%)
0.03
0.04
(d)
29
(b)
12
0.01
1. Provespa barthelemyi 2. Provespa anomala 3. Provespa nocturna 4. Vespula vulgaris 5. Vespula germanica 6. Dolichovespula media 7. Vespa basalis 8. Vespa binghami 9. Vespa soror 10. Vespa mandarinia 11. Vespa philippinensis 12. Vespa tropica 13. Vespa ducalis 14. Vespa fumida 15. Vespa crabro 16. Vespa dybowskii 17. Vespa mocsaryana 18. Vespa affinis 19. Vespa orientalis 20. Vespa analis 21. Vespa luctuosa 22. Vespa fervida 23. Vespa simillima 24. Vespa bicolor 25. Vespa bellicosa 26. Vespa multimaculata 27. Vespa vivax 28. Vespa velutina (Java) 29. V. velutina (Lombok) 30. V. velutina (Flores) 31. V. velutina (Sulawesi) 32. V. velutina (Thailand) 33. V. velutina (Vietnam) 34. V. velutina (Nepal) 35. V. velutina (Yunnan) 36. V. velutina (Zhejiang)
PC2 (19.07%) −0.04 −0.03 −0.02 −0.01 0.00
(a)
0.01
A. PERRARD ET AL.
PC2 (14.05%) −0.01 0.00
4
33 31 30
36 35 34 32
28
4 5 −0.02 0.00
−0.02
0.02 0.04 0.06 PC1 (54.7 %)
0.08
0.10
−0.01 0.00 PC1 (27.26%)
0.01
Fig. 2 Relationships between Vespinae species and populations and variation of the wing shapes in the multivariate tangent shape space. (a) Phylogeny of the Vespinae. (b–d) Two first principal components (PCs) of the variation of all Vespinae (b), Vespa (c) or Vespa velutina (d) wing shapes in the tangent shape space with the phylogeny mapped. Most of the groups appear overlapping because of the high dimensionality of the wing shape variation. Stars indicate nodes of clades with highly diverging wing shapes. Percentages indicate the amount of total shape variation expressed by the PC axes. Ancestral states of Vespa (c) and of V. velutina (d) are indicated by an arrow, and their relationships are depicted in dark grey.
Relationships between taxa were determined following the phylogenies of the genera Vespa and Provespa (Saito & Kojima, 2011; Perrard et al., 2013) combined with the genetic analysis of V. velutina populations (Perrard et al., 2014; Fig. 2a). Two main genetic lineages were distinguished within V. velutina, but relationships within these lineages were not well supported and are considered here as unresolved. Species and population values used for the phylogenetic signal test were computed as the mean shape of each taxonomic group. Phy-
logenetic signal was tested with 10 000 permutations on both the raw and allometry-corrected data for the whole data set and the data set of Vespa specimens. Allometry, the influence of size on shapes, may have an influence on the phylogenetic signal because some species groups or genera have similar sizes, which may potentially induce a similarity of shapes within these groups. The evolutionary allometry (Klingenberg, 1996) has been estimated by a multivariate regression of the taxa mean wing shapes on the taxa mean wing sizes, defined as the
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
5
Results
whereas PCAPop, based on more specimens, accounted for only 6.03% (Table 2). The dimensionality of the variation in the tangent space increased progressively from the Vespinae genera subsample to the V. velutina populations’ subsample (Fig. 4). The first PC of PCAGen explained more of the total variation than the ones of PCASp and of PCAPop. In addition, fewer variables of PCAGen were required to describe 95% of the variation (Table 2). These results could not be explained by the smaller sample size and fewer of groups of PCAGen: the dimensionality increase remained clear even with subsamples of PCASp and of PCAPop having the same sample size and group number than PCAGen.
Variability of the wing shape in the tangent space
Gaussian mixture models
More than half of the variation in the Vespinae specimens was contained in one dimension. It described the variation between the rounded forewings of Provespa and the more elongated forewings of the other species (Figs 2b and 3). The wing shape variability decreased with the level of classification: variation in PCAGen reached 36.35% of the variation in the entire sample,
Gaussian mixture models applied on the entire data set resulted in 30 clusters highly congruent with species groups: 78.92% of the 1025 specimens were in the cluster of their own species, all specimens being in a cluster of their own genus. For 17 of the 28 species, every specimen of a same species was gathered in one cluster. For six species, including V. velutina, more than 15% of their specimens were included in a second cluster. These second groups were specific for two species. Vespa velutina specimens were grouped into six clusters, two of them being shared with Vespa vivax and Vespa bicolor specimens. In the subsample of genera, five Gaussian clusters corresponded to the four genera, Provespa specimens being grouped in two clusters. In the subsample of the 22 Vespa species, 18 clusters were found. Most of these clusters grouped individuals from one species, resulting in 85% of the 720 Vespa specimens grouped with individuals of their own species (Fig. 5). When applied to the 338 specimens of V. velutina, GMM resulted in only one cluster of all specimens.
log-transformed centroid size of the wing. Phylogenetic signal for corrected allometric data was computed on the residuals of this regression. Every analysis was conducted with the software R (R Development Core Team, 2014) using the packages ‘Rmorph’ (Baylac, 2013), ‘mclust’ (Fraley & Raftery, 2002), ‘MASS’ (Venables & Ripley, 2002) and ‘ape’ (Paradis et al., 2004) as well as original R functions for the test of the phylogenetic signal. These functions were based on the algorithm described by McArdle & Rodrigo (1994) for ancestral states reconstruction of continuous characters minimizing squared parsimony.
(a)
(b)
Differences between the taxonomic groups
(c)
Fig. 3 Wing shape variation between Vespa and Provespa. (a) Difference in mean shapes of Vespa (grey) and Provespa (black) wings. (b) Wing of Vespa luctuosa. (c) Wing of Provespa barthelemyi. Scale bar is 1 mm.
The CVA with cross-validation on the whole sample returned 88.29% of correct assignment of specimens to their species and populations. It confirmed that 100% of the 1025 specimens could be attributed to the right genus and 97.07% to the right species. Only 70.71% of V. velutina specimens were reassigned to the right population, but the result is quite unbalanced between continental and insular populations: 96.12% (125 of 129) of the insular specimens were attributed to their correct island, whereas only 55.02% (115 of 209) of the continental specimens were attributed to their own population. Correct assignment tended to increase with the number of variables and with the exclusion of small groups, but the trends between taxonomic levels remained similar (see Supporting Information). Reassignment based on the k-nearest neighbours also
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
6
A. PERRARD ET AL.
Table 2 Wing shape variation in the multivariate tangent space for the entire Vespinae sample and for each taxonomic level subsample.
Total sample Vespinae genera Vespa species Vespa velutina populations
N indiv.
N groups
r2 (10 5)
N PC (95%)
1st PC (%)
2nd PC (%)
3rd PC (%)
1025 159 720 338
36 4 22 9
126.0 45.8 40.8 7.6
12 6 16 18
54.68 69.94 30.35 27.26
12.59 13.46 19.07 14.05
7.02 6.36 11.46 9.49
40
60
80
PVespa = 0) and allometry-corrected data (D2total = 9.133 9 10 3, Ptotal = 0; D2Vespa = 4.071 9 10 3, PVespa = 0). The phylogeny mapped in the tangent space showed that wing shapes tended to be grouped according to their clade in most cases. In the subsample of V. velutina, the first axis distinguished populations of the two main genetic lineages. Nonetheless, the influence of the phylogeny on the variation in the wing venation was not complete, as shown by the many overlapping lineages and the diverging nodes in the tangent space. As such, strong divergences were found between the wing shapes of Provespa and those of the other vespine genera, between the Vespa mandarinia species group and the Vespa tropica species group, as well as between the three closely related species Vespa affinis, Vespa orientalis and Vespa mocsaryana (Fig. 2b,c).
Vespinae genera Vespa species V. velutina populations
20
Percentage of variation
100 95
The column ‘N indiv.’ indicates the sample size, ‘N groups’ refers to the number of taxonomic groups compared, ‘r2’ refers to the variance in the total shape space, and ‘N PC’ refers to the number of principal components (or independent dimensions) required to explain 95% of the shape variation in the tangent space comprising 34 dimensions. The last three columns describe the percentage of variation explained by the three-first axes of the principal component analyses for each sample. The variance in the total shape space was assessed by the sum of the PC scores’ variances normalized by the ratio between each subsample size and the total sample size.
0
5
10 15 20 25 Number of dimensions
30
34
Fig. 4 Cumulative percentage of the total variation explained by the first components of principal component analyses on the Vespinae genera, Vespa species and Vespa velutina populations’ subsamples. Error bars are 95% estimates from 1000 bootstrap replicates.
confirmed these CVA results. The variation in wing shape means was 6.43 P times higher between genera than P between species ( 2genera = 0.00274, Ngenera = 4 2 and species = 0.00045, Nspecies = 22). The variation between species was 4.54 times higher than P that found between populations of V. velutina ( 2populations = 0.00009, Npopulations = 9). However, this variation between populations of the same species was slightly lower than the overall intrapopulation variation estimated on population-centred shape data in V. velutina P ( 2intrapopulation = 0.00015, Nintrapopulation = 338). Phylogenetic signal Phylogenetic signal was significant both in the whole data set and in the Vespa subsample, with both raw (D2total = 9.144 9 10 3, Ptotal = 0; D2Vespa = 3.902 9 10 3,
Discussion Wing venation shape as a taxonomic marker Our results showed that the evolution of the wing venation may explain its discriminating power by inducing structured, original and progressive differences of the venation shape between the different taxa. The genetic drift between taxa may have thus played a major role in the wing evolution. First, the main structure of the wing shape variation was found congruent with the taxonomic groups of social wasps. Many components could have been involved in the variation of wing shapes of 1025 specimens originating from a wide range of localities and belonging to species with largely overlapping distributions (Carpenter & Kojima, 1997). Gaussian mixture models showed, nonetheless, that the distribution of wing shapes in the tangent shape-space reflected species or groups of closely related species. These results confirmed that, within sex and caste, the genetic difference between species had a stronger influence on the wing shape than other biotic or abiotic factors. Secondly, the wing shape variation was highly multivariate: its evolution did not follow simple variation patterns even if the variation was not high in all 34
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Fig. 5 Species composition of the 18 Gaussian clusters based on the wing shape distribution of the Vespa subsample. Most clusters group mainly specimens from one species or from closely related species. The largest cluster (N = 77) groups many specimens from various species without taxonomical tendency. Species legends are as follows: a. Vespa bicolor; b. Vespa crabro; c. Vespa fervida; d. Vespa luctuosa; e. Vespa mandarinia; f. Vespa mocsaryana; g. Vespa multimaculata; h. Vespa philippinensis; i. Vespa simillima; j. Vespa soror; k. Vespa vivax; l. Vespa analis; m. Vespa ducalis; n. Vespa bellicosa.
Vespa affinis V. analis V. basalis V. analis V. bellicosa V. bicolor V. simillima V. binghami V. crabro V. ducalis V. dybowskii V. fervida V. fumida V. mandarinia V. mocsaryana V. multimaculata V. orientalis V. tropica V. soror V. velutina 0.0
0.2
dimensions of the tangent space (Fig. 4). This high dimensionality of the shape variation was previously reported in other studies on flies and honey bees (Monteiro et al., 2002; Mezey & Houle, 2005). Under a model of random evolution, high dimensionality of variation increases the potential ways of divergence during the evolutionary history. As a result, once two groups have diverged, they are less likely to converge to an exact same location again in the multidimensional space (C.P. Klingenberg, personal communication). The relatively high dimensionality of the wing shape variation may have promoted the evolution of specific venation shapes by preventing complete convergence of the wing between different taxa. Furthermore, the wing shape variation not only reflected the taxonomic groups, but was also partially congruent with the evolutionary history shared by these groups. The CVA showed that the distribution of wing shapes of a single genus in the tangent space was not overlapping with the ones from its closest genera. On the other hand, the distinction was lower between species: 2.93% of the individuals were attributed to a wrong species within the same genus and even lower between populations of the single species V. velutina. The insular populations of V. velutina were more readily discriminated than the continental populations (96.12% vs. 55.02%). It is therefore possible that genetic isolation experienced in islands enhanced the evolution of distinct wing shapes (Camara et al., 2006). Furthermore, both the GMM and CVA misattributions between species mainly involved closely related species, suggesting that closely related species tended to have more similar wing shapes. The comparison of variances at the different levels of the classification confirmed that the variation between groups increased with the taxonomic level. Finally, the significance of the phylogenetic signal showed that the variation in wing shape
0.4
c. a. b.
7
n = 34 V. crabro n = 19 g. k. n = 27 V. velutina n = 77 f. h. i. j. V. velutina V. vivax n = 64 n = 24 n = 47 V. dybowskii n = 30 n = 20 n = 43 V. bellicosa n = 13 n = 22 n = 26 V. bicolor V. bellicosa l. n = 47 n = 38 V. philippinensis m. n. n = 43 n = 32 V. mandarinia n = 22
d. e.
V. luctuosa
0.6
0.8
1.0
between the populations, species and genera was congruent with their evolutionary history. The unresolved nodes in the phylogeny did not interfere with this result: by increasing the length of the tree in the tangent space (Klingenberg & Gidaszewski, 2010), unresolved nodes decrease the probability of a significant signal. The phylogenetic signal being based on a Brownian motion model of the phenotype evolution, this last result supports an evolution of the phenotype mainly driven by genetic drift. Altogether, these results are further evidence that the variation in the wing venation is largely shaped by the evolutionary history and could be used to identify the genera and species of insects presenting a venation pattern suitable for landmark analyses. Other factors of wing shape variation Our results also suggested that drift may not have been the only factor influencing the wing shape evolution in social wasps. Although there is a significant phylogenetic signal in the overall variation, some of our results were more likely explained by the influence of external factors than by pure drift.
Body shape The main wing shape variation found in our sampling was the difference between the apically rounded Provespa wings and the more elongated wings of other Vespinae (Figs 2b and 3). A similar variation of the wing was found in Polistinae where rounded wings were related to small body sizes (Garcıa & Sarmiento, 2012). In our case, this relationship does not hold as Provespa specimens were larger than those of Vespula and Dolichovespula. Some of the main differences between Provespa and other Vespinae genera such as enlarged ocelli and brown coloration are adaptations to nocturnal hab-
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
A. PERRARD ET AL.
(a)
0.01
1 v4 v1 2
v2 v3 v5
PC2 −0.01
0.00
3
v7 6
v6 4 v8 5
v9 7
−0.02
its. However, no straightforward link can be established between wing shape and nocturnal habits. Furthermore, one hornet species, Vespa binghami, is nocturnal and shares the enlarged ocelli and brown coloration with Provespa while retaining a Vespa-like wing shape. The divergent wings of Provespa species may thus be more likely related to their slender and proportionally lighter body than those of yellow jackets and hornets. The wing venation has a role in wing rigidity (Combes & Daniel, 2003) which is important for flight (Young et al., 2009; Mountcastle & Combes, 2013). Furthermore, the shape of the wing is likely involved in the production of the aerodynamic forces during the flight (Young et al., 2009). A variation in the ratio of the body mass compared with its length during the evolution of Provespa wasps may have induced a shift in the optimal wing shape. These divergent selective pressures could have enhanced the development of divergent wing shapes between Provespa species and the other Vespinae.
8
−0.04
8
−0.02
−0.01
0.01 0.00 PC1
0.02
0.03
(b) PC2
Climate The wing shape of V. orientalis, the only species of hornet living in arid climates, was distinguished from those of all other Vespa species by the second PC of variation within the genus (Fig. 6). Furthermore, the same axis of variation separated the mean shapes of continental and Indonesian populations of V. velutina, the latter being closer to the mean shapes of the four Indonesian species of the same clade than V. velutina. On the other hand, continental populations of V. velutina were closer to continental species of the clade including V. velutina (Fig. 6a). These results suggested that the second PC of the variation within Vespa described a shape change influenced by geographical and potentially climatic components. Latitudinal and altitudinal clines of wing variation were found in several species of Drosophila (e.g. Dahlgaard et al., 2001; Sambucetti et al., 2006; Pitchers et al., 2013). Furthermore, wing shape development is known to be influenced by temperature (Debat et al., 2006; Pitchers et al., 2013). The variation of climatic parameters such as humidity and temperature in different geographical ranges may have influenced the shape variation between wasp species either through plasticity of the phenotype or by producing different selective pressures. Phenotypic plasticity, the ability to develop different phenotypes with a same genetic background in response to different environments, has been previously found in wings of Drosophila in relation to thermal clines (Debat et al., 2006; Pitchers et al., 2013). Phenotypic plasticity is likely to have played a role in the observed variation of Vespa wings because this variation was detected between species as well as within a single species (Fig. 3). However, hornets develop in enclosed nests with buffered temperature which should limit the effect of climatic differences on wings develop-
Fig. 6 Variation in wing shape potentially influenced by climatic factors. (a) Two first dimensions of the tangent shape space of Vespa wings with the phylogeny mapped. The second principal component (PC) differentiates continental (white squares) and Indonesian (white circles) species and populations of the clade of Vespa velutina as well as the only species from arid zones Vespa orientalis (star). 1: Vespa vivax; 2: Vespa simillima; 3: Vespa bicolor; 4: Vespa multimaculata; 5: Vespa bellicosa; 6: Vespa luctuosa; 7: Vespa fervida; 8: V. orientalis; v1–9: V. velutina. v1: Vietnam; v2: Zhejiang; v3: Yunnan; v4: Thailand; v5: Nepal; v6: Java; v7: Sulawesi; v8: Flores; v9: Lombok. (b) Wing shape variation along the second PC axis from negative (grey) to positive (black) values.
ment (Jones & Oldroyd, 2006). Furthermore, the geographical clines could not be entirely explained by external factors even in Drosophila (Dahlgaard et al., 2001; Sambucetti et al., 2006). Part of the wing variation related to the geography in hornets, such as the highly divergent wings of V. orientalis, may thus come from fixed genetic differences driven by divergent selective pressures. Phenotypic plasticity promotes evolutionary answers to such selective pressures through ‘genetic accommodation’ (West-Eberhard, 2005). The combination of phenotypic plasticity and divergent selection may thus explain the observed variation of wing shape congruent with the geography.
Mimicry The wing shape of V. affinis, a species closely related to V. mocsaryana, had an extreme value on the first PC of the PCASp, close to the shape of its distant relative V. tropica (Fig. 2c). These two species present also similar geo-
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
graphical distribution and coloration, a dark body with a pale area on the two-first metasomal terga. On the other hand, wing shapes of V. mocsaryana and Vespa basalis were the closest, both species having reddish colour and almost superimposed distributional range (Carpenter & Kojima, 1997). The hypothesis of a convergence due to shared climatic conditions cannot be retained for these species as these similar distributions, and thus climatic areas, are largely overlapping with those of the other Vespa species. However, wing shape was related to wing sound production in Drosophila (Routtu et al., 2007; Menezes et al., 2013), and buzzing was proven to enhance the aposematism efficiency in wasps (Hauglund et al., 2006). Wing shape could thus be selected by mimicry pressures on the flight sound. A second explanation could be that genes influencing the patterning of pigmentation, such as genes Hox (Kopp et al., 1999; Perrard et al., 2014), may also be involved in the development of similar venation patterns. Further research is required to properly test the potential influence of these factors on wing shape and to understand the origin of the resulting variation, from plasticity or fixed genetic differences. Our results, nonetheless, demonstrate that these external influences are not strong enough to blur the wing shape differences due to the divergence between genera and species.
Conclusion Wing shape measurement and analysis is fast and requires only a reference data set, a camera and a computer. Previous studies have shown that this is an efficient identification criterion on a wide variety of insects. These include economically important species such as vectors of diseases (Jirakanjanakit & Dujardin, 2005; Henry et al., 2010), pests (Schutze et al., 2011) and pollinators (Francoy et al., 2008, 2009). Our results provide new insights confirming and explaining the efficiency of the venation shape for insect identification: the main source of variation in the venation shape was the divergence between taxa, and this divergence occurred progressively, with rare and incomplete secondary convergence. With an initial effort to produce the reference data sets, using identified specimens from natural history collections, for example, it could become a quick and inexpensive identification tool (Steinhage et al., 2001; Camara et al., 2006; Schutze et al., 2011). Furthermore, standardized wing imaging can be performed using portable material such as a camera with a tripod and microscope blades (Perrard et al., 2012). Wing shape could thus enable insect identification directly on the field, which would represent an important advance for many ecological studies and pest surveys. Finally, our results not only explained the efficiency of the wing shape as an identification criterion, but also showed clues to other factors influencing its evolution.
9
The wing venation has a clear functional role in flight by determining the rigidity of the wing (Combes & Daniel, 2003; Young et al., 2009; Mountcastle & Combes, 2013). However, the functional importance of the venation shape is still poorly understood (Routtu et al., 2007; Mountcastle & Combes, 2013). In our analysis, the wing shape variation tended to be highly multivariate and homogeneous within a species while shapes diverged along specific directions between distant species or genera. Whether these evolutionary pathways in the space of possible shapes resulted from developmental constraints or selective pressures is yet to be elucidated. Further analyses of these patterns of shape variation in regard to potential adaptive processes could help to understand the role of the wing shape and venation in insect flight.
Acknowledgments This work was funded by two ATM (‘Formes possibles, formes r ealis ees’ and ‘Biodiversit e actuelle et fossile’) grants from the Museum National d’Histoire Naturelle, the Bourse Eole of the Franco-Netherland Alliance and the Bourse Germaine Cousin of the Soci et e Entomologique de France. We thank C.P. Klingenberg, P. Mitteroecker and the anonymous reviewer for their helpful comments on the manuscript. We are grateful to J. Haxaire, J.-L. Renesson and P. Tripotin who shared their private collections. Special thanks are due to curators A. Drumont and G. Yvonnick (Institut Royal des Sciences Naturelles de Belgique, Brussels), G. Broad (Natural History Museum, London), F. and J. Gusenleitner (Ober€ osterreiches Landesmuseum, Linz), K. van Achterberg (Nationaal Natuurhistorisch Museum Naturalis, Leiden), H. Nugroho (Museum Zoologicum Bogoriense, Bogor) and X. Chen (University of Zhejiang, Hangzhou) who allowed to study in their respective museum collections. This study also benefited from the public collections of the American Museum of Natural History (New York, USA) and the Mus eum National d’Histoire Naturelle (MNHN, Paris, France). This study is a contribution of the Morphometric facility of the MNHN (plate-forme Morphom etrie, CNRS-MNHN UMS 2700).
Data accessibility Landmark coordinates of the measured wings, R scripts and R functions to compute the phylogenetic signals are available from the Dryad Digital Repository: http://doi. org/10.5061/dryad.2cv1h.
References Adams, D.C., Rohlf, F.J. & Slice, D.E. 2004. Geometric morphometrics: ten years of progress following the ‘revolution’. Ital. J. Zool. 71: 5–16.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
10
A. PERRARD ET AL.
Bai, M., McCullough, E., Song, K.-Q., Liu, W.-G. & Yang, X.-K. 2011. Evolutionary constraints in hind wing shape in Chinese dung beetles (Coleoptera: Scarabaeinae). PLoS ONE 6: e21600. Baylac, M. 2013. Rmorph: a R geometric and multivariate morphometrics library. Available from the author: [email protected] Baylac, M. & Daufresne, T. 1996. Wing venation variability in Monarthropalpus buxi (Diptera, Cecidomyiidae) and the quaternary coevolution of box (Buxus sempervirens L.) and its midge: a geometrical morphometric analysis. In: Advances in Morphometrics (L. Marcus, ed.), pp. 285–302. Plenum Press, New York, NY. Baylac, M., Villemant, C. & Simbolotti, G. 2003. Combining geometric morphometrics with pattern recognition for the investigation of species complexes. Biol. J. Linn. Soc. 80: 89– 98. Bickford, D., Lohman, D.J., Sodhi, N.S., Ng, P.K., Meier, R., Winker, K. et al. 2007. Cryptic species as a window on diversity and conservation. Trends Ecol. Evol. 22: 148–155. Bookstein, F.L. 1991. Morphometric Tools for Landmark Data: Geometry and Biology. Cambridge University Press, New York, NY. Camara, M., Caro-Ria~ no, H., Ravel, S., Dujardin, J., Hervouet, J., De Mee€ us, T. et al. 2006. Genetic and morphometric evidence for population isolation of Glossina palpalis gambiensis (Diptera: Glossinidae) on the Loos islands, Guinea. J. Med. Entomol. 43: 853–860. Carpenter, J.M. & Kojima, J. 1997. Checklist of the species in the subfamily Vespinae (Insecta: Hymenoptera: Vespidae). Nat. Hist. Bull. Ibaraki Univ. 1: 51–92. Combes, S.A. & Daniel, T.L. 2003. Flexural stiffness in insect wings I. Scaling and the influence of wing venation. J. Exp. Biol. 206: 2979–2987. Comstock, J.H. & Needham, J.G. 1899. The wings of insects. Am. Nat. 33: 845–860. Dahlgaard, J., Hasson, E. & Loeschcke, V. 2001. Behavioral differentiation in oviposition activity in Drosophila buzzatii from highland and lowland populations in Argentina: plasticity or thermal adaptation? Evolution 55: 738–747. De Meulemeester, T., Michez, D., Aytekin, A.M. & Danforth, B.N. 2012. Taxonomic affinity of halictid bee fossils (Hymenoptera: Anthophila) based on geometric morphometrics analyses of wing shape. J. Syst. Palaeontol. 10: 755–764. Debat, V., Milton, C.C., Rutherford, S., Klingenberg, C.P. & Hoffmann, A.A. 2006. Hsp90 and the quantitative variation of wing shape in Drosophila melanogaster. Evolution 60: 2529– 2538. Dryden, I.L. & Mardia, K.V. 1998. Statistical Shape Analysis. John Wiley & Sons, Chichester. Fraley, C. & Raftery, A.E. 1998. How many clusters? Which clustering method? Answers via model-based cluster analysis. Comput. J. 41: 578–588. Fraley, C. & Raftery, A.E. 2002. Model-based clustering, discriminant analysis, and density estimation. J. Am. Stat. Assoc. 97: 611–631. Francoy, T.M., Wittmann, D., Drauschke, M., M€ uller, S., Steinhage, V., Bezerra-Laure, M.A. et al. 2008. Identification of Africanized honey bees through wing morphometrics: two fast and efficient procedures. Apidologie 39: 488–494. Francoy, T.M., Silva, R.A.O., Nunes-Silva, P., Menezes, C. & Imperatriz-Fonseca, V.L. 2009. Gender identification of five
genera of stingless bees (Apidae, Meliponini) based on wing morphology. Genet. Mol. Res. 8: 207–214. Garcıa, Z. & Sarmiento, C.E. 2012. Relationship between body size and flying-related structures in neotropical social wasps (Polistinae, Vespidae, Hymenoptera). Zoomorphology 131: 1– 11. Gaston, K.J. & O’Neill, M.A. 2004. Automated species identification: why not? Philos. Trans. R. Soc. Lond. B Biol. Sci. 359: 655–667. Grimaldi, D.A. & Engel, M.S. 2005. Evolution of the Insects. Cambridge University Press, Cambridge. Hauglund, K., Hagen, S.B. & Lampe, H.M. 2006. Responses of domestic chicks (Gallus gallus domesticus) to multimodal aposematic signals. Behav. Ecol. 17: 392–398. Henry, A., Thongsripong, P., Fonseca-Gonzalez, I., Jaramillo-Ocampo, N. & Dujardin, J.-P. 2010. Wing shape of dengue vectors from around the world. Infect. Genet. Evol. 10: 207–214. Huey, R.B. & Bennett, A.F. 1987. Phylogenetic studies of coadaptation: preferred temperatures versus optimal performance temperatures of lizards. Evolution 41: 1098–1115. Jirakanjanakit, N. & Dujardin, J.-P. 2005. Discrimination of Aedes aegypti (Diptera: Culicidae) laboratory lines based on wing geometry. SE. Asian J. Trop. Med. 36: 858–861. Jones, J.C. & Oldroyd, B.P. 2006. Nest thermoregulation in social insects. Adv. In Insect Phys. 33: 153–191. Kendall, D.G. 1977. The diffusion of shape. Adv. Appl. Probab. 9: 428–430. Klingenberg, C.P. 1996. Multivariate allometry. In: Advances in Morphometrics (L. Marcus, ed.), pp. 23–49. Plenum Press, New York, NY. Klingenberg, C.P. & Gidaszewski, N.A. 2010. Testing and quantifying phylogenetic signals and homoplasy in morphometric data. Syst. Biol. 59: 245–261. Klingenberg, C.P., Barluenga, M. & Meyer, A. 2002. Shape analysis of symmetric structures: quantifying variation among individuals and asymmetry. Evolution 56: 1909–1920. Kopp, A., Blackman, R.K. & Duncan, I. 1999. Wingless, decapentaplegic and EGF receptor signaling pathways interact to specify dorso-ventral pattern in the adult abdomen of Drosophila. Development 126: 3495–3507. Marcus, L., Hingst-Zaher, E. & Zaher, H. 2000. Application of landmark morphometrics to skulls representing the orders of living mammals. Hystrix 11: 27–47. Matias, A., de la Riva, J., Torrez, M. & Dujardin, J.-P. 2001. Rhodnius robustus in Bolivia identified by its wings. Mem. I. Oswaldo Cruz 96: 947–950. Mayr, E. 1969. Principles of Systematic Zoology. McGraw-Hill, New York, NY. McArdle, B. & Rodrigo, A.G. 1994. Estimating the ancestral states of a continuous-valued character using squared-change parsimony: an analytical solution. Syst. Biol. 43: 573–578. Meier, R., Shiyang, K., Vaidya, G. & Ng, P.K. 2006. DNA barcoding and taxonomy in Diptera: a tale of high intraspecific variability and low identification success. Syst. Biol. 55: 715–728. Menezes, B.F., Vigoder, F.M., Peixoto, A.A., Varaldi, J. & Bitner-Mathe, B.C. 2013. The influence of male wing shape on mating success in Drosophila melanogaster. Anim. Behav. 85: 1217–1223. Mezey, J.G. & Houle, D. 2005. The dimensionality of genetic variation for wing shape in Drosophila melanogaster. Evolution 59: 1027–1038.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Mitteroecker, P. & Bookstein, F. 2011. Linear discrimination, ordination, and the visualization of selection gradients in modern morphometrics. Evol. Biol. 38: 100–114. Monteiro, L.R., Diniz-Filho, J.A.F., Reis, S.F. & Ara ujo, E.D. 2002. Geometric estimates of heritability in biological shape. Evolution 56: 563–572. Mountcastle, A.M. & Combes, S.A. 2013. Wing flexibility enhances load-lifting capacity in bumblebees. Proc. Roy. Soc. B: Biol. Sci. 280: 20130531. Outomuro, D., Adams, D.C. & Johansson, F. 2013. The evolution of wing shape in ornamented-winged damselflies (Calopterygidae, Odonata). Evol. Biol. 40: 300–309. Paradis, E., Claude, J. & Strimmer, K. 2004. APE: analyses of phylogenetics and evolution in R language. Bioinformatics 20: 289–290. Perrard, A., Villemant, C., Carpenter, J.M. & Baylac, M. 2012. Differences in caste dimorphism among three hornet species (Hymenoptera: Vespidae): forewing size, shape and allometry. J. Evol. Biol. 25: 1389–1398. Perrard, A., Pickett, K.M., Villemant, C., Kojima, J. & Carpenter, J.M. 2013. Phylogeny of hornets: a total evidence approach (Hymenoptera, Vespidae, Vespinae, Vespa). J. Hymenopt. Res. 32: 1–15. Perrard, A., Arca, M., Rome, Q., Muller, F., Tan, J., Bista, S. et al. 2014. Geographic variation of melanisation patterns in a hornet species: genetic differences, climatic pressures or aposematic constraints? PLoS ONE 9: e94162. Pitchers, W., Pool, J.E. & Dworkin, I. 2013. Altitudinal clinal variation in wing size and shape in African Drosophila melanogaster: one cline or many? Evolution 67: 438–452. R Development Core Team. 2014. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. http://R-project.org (accessed 16 March 2014). Ripley, B.D. 1996. Pattern Recognition and Neural Networks. Cambridge University Press, Cambridge. Rohlf, F.J. 2010. tpsDig v2. 16. Department of Ecology and Evolution, State University of New York, Stony Brook, NY (Distributed by author). Routtu, J., Mazzi, D., Van Der Linde, K., Mirol, P., Butlin, R.K. & Hoikkala, A. 2007. The extent of variation in male song, wing and genital characters among allopatric Drosophila montana populations. J. Evol. Biol. 20: 1591–1601. Saito, F. & Kojima, J. 2011. Phylogenetic analysis and biogeography of the nocturnal hornets, Provespa (Insecta: Hymenoptera: Vespidae: Vespinae). Species Div. 16: 65–74. Sambucetti, P., Loeschcke, V. & Norry, F.M. 2006. Developmental time and size-related traits in Drosophila buzzatii along an altitudinal gradient from Argentina. Hereditas 143: 77–83.
11
Schutze, M.K., Jessup, A. & Clarke, A.R. 2011. Wing shape as a potential discriminator of morphologically similar pest taxa within the Bactrocera dorsalis species complex (Diptera: Tephritidae). Bull. Entomol. Res. 102: 103–111. Steinhage, V., Arbuckle, T., Schr€ oder, S., Cremers, A.B. & Wittmann, D. 2001. ABIS: automated identification of bee species. In: BIOLOG. Workshop, German Programme on Biodiversity and Global Change Status Report (W.G. Berendsohn, ed.), pp. 194–195. Bonn. Venables, W.N. & Ripley, B.D. 2002. Modern Applied Statistics with S, 4th edn. Springer Verlag, Berlin. Villemant, C., Simbolotti, G. & Kenis, M. 2007. Discrimination of Eubazus (Hymenoptera, Braconidae) sibling species using geometric morphometrics analysis of wing venation. Syst. Entomol. 32: 625–634. West-Eberhard, M.J. 2005. Developmental plasticity and the origin of species differences. Proc. Natl. Acad. Sci. USA 102: 6543–6549. Young, J., Walker, S.M., Bomphrey, R.J., Taylor, G.K. & Thomas, A.L.R. 2009. Details of insect wing design and deformation enhance aerodynamic function and flight efficiency. Science 325: 1549–1552. Zelditch, M., Swiderski, D.L., Sheets, H.D. & Fink, W.L. 2004. Geometric Morphometrics for Biologists: A Primer. Elsevier Academic Press, London.
Supporting information Additional Supporting Information may be found in the online version of this article: Table S1 Classification to genera. Table S2 Classification to species. Table S3 Classification to species excluding groups of sample sizes < 10. Table S4 Classification to species excluding groups of sample sizes < 15. Table S5 Classification to species excluding groups of sample sizes < 20. Table S6 Classification to populations. Table S7 Classification to populations excluding groups of sample size < 10. Table S8 Classification to populations excluding groups of sample size < 34. Received 29 May 2014; revised 18 September 2014; accepted 26 September 2014
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Evolution of wing shape in hornets: why is the wing venation efficient for species identification? A. PERRARD*†, M. BAYLAC†, J. M. CARPENTER* & C. VILLEMANT† *Division of Invertebrate Zoology, American Museum of Natural History, New York, NY, USA †UMR 7205 ISYEB, Museum National d’Histoire Naturelle, Paris, France
Keywords:
Abstract
Gaussian mixture models; genus Vespa; insect wing; phylogenetic signal; species recognition; tangent shape space.
Wing venation has long been used for insect identification. Lately, the characterization of venation shape using geometric morphometrics has further improved the potential of using the wing for insect identification. However, external factors inducing variation in wing shape could obscure specific differences, preventing accurate discrimination of species in heterogeneous samples. Here, we show that interspecific difference is the main source of wing shape variation within social wasps. We found that a naive clustering of wing shape data from taxonomically and geographically heterogeneous samples of workers returned groups congruent with species. We also confirmed that individuals can be reliably attributed to their genus, species and populations on the basis of their wing shape. Our results suggested that the shape variation reflects the evolutionary history with a potential influence of other factors such as body shape, climate and mimicry selective pressures. However, the high dimensionality of wing shape variation may have prevented absolute convergences between the different species. Wing venation shape is thus a taxonomically relevant marker combining the accuracy of quantitative characters with the specificity required for identification criteria. This marker may also highlight adaptive processes that could help understand the wing’s influence on insect flight.
Introduction Species identification and classification is a first and crucial step in many biological studies. It allows researchers to attribute information to recognizable entities (Mayr, 1969; Gaston & O’Neill, 2004). However, species recognition is often challenged by taxonomic diversity and the subtlety of interspecific variation (Meier et al., 2006; Bickford et al., 2007). With hundreds of thousands of distinct species, insects are one of the most diverse taxa of the living world (Grimaldi & Engel, 2005). Wing venation patterns are classically used for their identification at the order or family levels (Comstock & Needham, 1899). More recently, the quantification of the wing shape through geometric morphometrics became widely Correspondence: Adrien Perrard, Division of Invertebrate Zoology, American Museum of Natural History, 79th Street at Central Park West, New York, NY 10024, USA. Tel.: +1 212 496 3446; fax: +1 212 769 5277; e-mail: [email protected]
used to study subtle variation among insect genera, species and even populations, sex and castes within the same species (e.g. Baylac & Daufresne, 1996; Baylac et al., 2003; Villemant et al., 2007; Bai et al., 2011; De Meulemeester et al., 2012; Perrard et al., 2012; Outomuro et al., 2013). These previous studies showed that wing shape differs between most taxonomic groups and can be used for insect identification (Matias et al., 2001; Steinhage et al., 2001; Francoy et al., 2008, 2009; Henry et al., 2010). However, they did not address how the wing shape varies along a lineage: Is this marker pertinent for a generalized identification, independent from the taxonomic scale of the sampling? Is the morphological divergence between taxa the main source of variation in the venation shape? Are there other factors influencing its evolution? Geometric morphometrics aim at quantifying geometric shapes, a spatial description of objects excluding size, orientation and position (Bookstein, 1991; Adams et al., 2004). Insect wing shapes are usually described by
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
1
2
A. PERRARD ET AL.
landmarks placed at the vein intersections (Baylac & Daufresne, 1996; Klingenberg et al., 2002). Superimposition methods enable comparison of different shapes by extracting a standardized shape description from the raw landmark coordinates. The space of all possible shapes described by the same number of landmarks is a nonlinear multidimensional manifold called the ‘shape space’ (Kendall, 1977). Because most multivariate statistical analyses are designed for linear spaces, the aligned coordinates resulting from the superimposition are usually projected onto a linear tangent space. In biology, the difference between the shape space and the tangent space is considered negligible: biological variation in shape being very limited compared to the geometric possible shapes, the subspace of biological shapes remains close to the tangent space (Marcus et al., 2000; Zelditch et al., 2004). An ensemble of shapes can thus be represented as a cloud of points in the multivariate tangent space. To be a generalized marker for taxon identification, the variation in the wing venation should enable identification at different levels of the classification. It should be highly structured by evolution and should have undergone low convergence between taxa: the wing venation shape should thus be more influenced by the genetic drift between taxa than by other factors of variation. This implies several predictions as to the distribution of wing shapes within the tangent space: first, wing shapes should be distinct between taxa at different taxonomic levels; second, for a given taxonomic level, the intragroup variation should be lower than the intergroup variation. Finally, as a corollary of the two previous statements, the shape difference between distantly related groups (e.g. different genera)
is expected to be higher than between more closely related taxa. The average shape difference between taxa should thus reflect their classification level and, indirectly, their degree of relatedness. To test these hypotheses, we analysed the evolution of the wing shape of social wasps of the subfamily Vespinae at different classification levels, from genera to populations. Wing shape is influenced by external factors affecting the development of organisms (Debat et al., 2006), and this variation could blur the difference between groups, especially in case of phenotypic convergence. To address the full variability of shapes within each group and to test the difference between groups of the same taxonomic level, we used wild specimens collected from different localities for each species or populations. We first explored the structure of the variability of wing shapes in the tangent space using Gaussian-based clustering. Then, we estimated the difference between the taxonomic groups using linear models and cross-validation attribution tests. Finally, we tested the phylogenetic signal of the wing shape to establish whether the variation in wing shape between the taxonomic groups reflected their relatedness.
Materials and methods Social wasp specimens were gathered from natural history collections of eight public institutions and three private collections. We sampled 1025 specimens among the four genera of the subfamily Vespinae with a special emphasis on hornets by including the 22 species of the genus Vespa and nine populations of the yellow-legged hornet Vespa velutina Lepeletier, 1836 (Table 1).
Species
N
Species [Population]
N
Dolichovespula media (Retzius, 1783) Provespa anomala (de Saussure, 1854) Provespa barthelemyi (du Buysson, 1905) Provespa nocturna v. d. Vecht, 1935 Vespa affinis (Linnaeus, 1764) Vespa analis Fabricius, 1775 Vespa basalis Smith, 1852 Vespa bellicosa de Saussure, 1854 Vespa bicolor Fabricius, 1787 Vespa binghami du Buysson, 1905 Vespa crabro Linnaeus, 1758 Vespa ducalis Smith, 1852 , 1884 Vespa dybowskii Andre Vespa fervida Smith, 1858 Vespa fumida v. d. Vecht, 1959 Vespa luctuosa de Saussure, 1854 Vespa mandarinia Smith, 1852 Vespa mocsaryana du Buysson, 1905
10 51 9 11 33 34 27 32 25 24 49 33 21 23 13 19 30 32
rez, 1910 Vespa multimaculata Pe Vespa orientalis Linnaeus, 1771 Vespa philippinensis de Saussure, 1854 Vespa simillima Smith, 1868 Vespa soror du Buysson, 1905 Vespa tropica (Linnaeus, 1758) Vespa velutina Lepeletier, 1836 [Java] V. velutina [Lombok] V. velutina [Flores] V. velutina [Nepal] V. velutina [Sulawesi] V. velutina [Thailand] V. velutina [Vietnam] V. velutina [Yunnan] V. velutina [Zhejiang] Vespa vivax Smith, 1870 Vespula germanica (Fabricius, 1793) Vespula vulgaris (Linnaeus, 1758)
34 38 10 30 34 31 54 40 25 68 10 33 10 35 63 6 15 13
Table 1 Sampling of the species and populations of Vespinae.
Column ‘N’ indicates the sample size.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Right forewings of specimens were measured using the protocol for pinned specimens described in Perrard et al. (2012). Only workers were included to avoid caste dimorphism. Forewing shape was assessed using 19 landmarks recorded with TPSDig2 software (Rohlf, 2010; Fig. 1). As the tangent space is defined by the superimposition, landmarks of all the specimens were superimposed only once using a generalized Procrustes analysis and thus projected in the same tangent space (Dryden & Mardia, 1998). In this study, wing shapes were described as the scores of the 34 principal components (PC) that have eigenvalues > 0 of a principal component analysis (PCA) of these projections. They are referred to as ‘specimens’ in the text. Variability of wing shape in the tangent space We first explored the variability of the wing shape in the tangent space at three classification levels, from population to genera. Because of the sampling ranging from V. velutina populations to Vespinae genera, the species V. velutina and the genus Vespa were overrepresented in the main PCA compared with the other species and genera. To avoid the bias towards these groups when analysing the variation at the different levels of classification, we performed three additional PCA. The first included the 338 specimens of V. velutina from the nine populations sampled (PCAPop). The second included 50 specimens of V. velutina randomly selected among the 338 specimens and every specimen of the 21 other species of Vespa (PCASp). The third PCA included 50 specimens of Vespa randomly selected from the previous sample and representatives of the three other Vespinae genera: Dolichovespula, Provespa and Vespula (PCAGen). We compared the variances of these subsamples and of the total sample using the respective sums of the 34 PC scores’ variances normalized by the ratio between each subsample size and the total sample size.
3
Gaussian mixture models We looked for groups of similar shapes in the whole sample and in the different subsamples using a Bayesian clustering method called Gaussian mixtures models (GMM; Fraley & Raftery, 1998). This method requires no a priori information on groups or number of groups; it clusters individuals based on the likelihood of groups of Gaussian distributions in the tangent space. If the divergence between taxonomic groups was the main cause of shape variation, Gaussian clusters should be congruent with these groups. Differences between taxa The difference in wing shapes between groups was estimated by a reassignment method based on canonical variate analyses (CVA) with leave-one-out cross-validation. Under the hypothesis of a complete separation of the groups within the tangent space, a perfect attribution would be expected with cross-validated reassignment. Because several groups had small sample sizes compared with the number of variables, the validity of CVA classifications was investigated (Mitteroecker & Bookstein, 2011). Reassignments were performed with different number of variables and excluding progressively the least sampled groups. The results were compared to reassignment obtained from the k-nearest neighbours nonparametric method with k = 1, 3 and 5 (Ripley, 1996). The amount of variation at the different classification levels was determined on the whole sample by the sum of variances (∑2) between the mean shapes of (i) the four Vespinae genera; (ii) the 22 Vespa species; and (iii) the nine populations of V. velutina. Under the hypothesis of a neutral evolution of the shape driven by genetic drift, the sum of squares of mean shapes of higher classification levels should be higher than those corresponding to lower levels or intragroup variation. Phylogenetic signal
Fig. 1 Location of the 2D landmarks measured on the right forewing of specimens.
The congruence between the wing shape variation and the evolutionary history was tested using the phylogenetic signal test from Klingenberg & Gidaszewski (2010). The phylogenetic tree was mapped within the tangent space by computing the ancestral states minimizing the squared parsimony (Huey & Bennett, 1987). The squared length of the tree within the tangent space (D2) was then used to test the phylogenetic signal with a significance level of 0.05. Under a Brownian motion model of evolutionary change of the wing, the phylogeny should be significantly shorter than random trees in the tangent space. The null distribution of the squared tree lengths was computed by permutation of the taxa across the phylogeny.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
c
d
0.02
b
(c)
d 17
14
7
13
27 23 8 24 9 20 21 16 26 10 25 15 22
11
0.02 PC2 (12.6 %) −0.04 −0.02 0.00
3 12 c 6
18
19
−0.02 −0.01 0.00 0.01 0.02 PC1 (30.35%)
0.03
0.04
(d)
29
(b)
12
0.01
1. Provespa barthelemyi 2. Provespa anomala 3. Provespa nocturna 4. Vespula vulgaris 5. Vespula germanica 6. Dolichovespula media 7. Vespa basalis 8. Vespa binghami 9. Vespa soror 10. Vespa mandarinia 11. Vespa philippinensis 12. Vespa tropica 13. Vespa ducalis 14. Vespa fumida 15. Vespa crabro 16. Vespa dybowskii 17. Vespa mocsaryana 18. Vespa affinis 19. Vespa orientalis 20. Vespa analis 21. Vespa luctuosa 22. Vespa fervida 23. Vespa simillima 24. Vespa bicolor 25. Vespa bellicosa 26. Vespa multimaculata 27. Vespa vivax 28. Vespa velutina (Java) 29. V. velutina (Lombok) 30. V. velutina (Flores) 31. V. velutina (Sulawesi) 32. V. velutina (Thailand) 33. V. velutina (Vietnam) 34. V. velutina (Nepal) 35. V. velutina (Yunnan) 36. V. velutina (Zhejiang)
PC2 (19.07%) −0.04 −0.03 −0.02 −0.01 0.00
(a)
0.01
A. PERRARD ET AL.
PC2 (14.05%) −0.01 0.00
4
33 31 30
36 35 34 32
28
4 5 −0.02 0.00
−0.02
0.02 0.04 0.06 PC1 (54.7 %)
0.08
0.10
−0.01 0.00 PC1 (27.26%)
0.01
Fig. 2 Relationships between Vespinae species and populations and variation of the wing shapes in the multivariate tangent shape space. (a) Phylogeny of the Vespinae. (b–d) Two first principal components (PCs) of the variation of all Vespinae (b), Vespa (c) or Vespa velutina (d) wing shapes in the tangent shape space with the phylogeny mapped. Most of the groups appear overlapping because of the high dimensionality of the wing shape variation. Stars indicate nodes of clades with highly diverging wing shapes. Percentages indicate the amount of total shape variation expressed by the PC axes. Ancestral states of Vespa (c) and of V. velutina (d) are indicated by an arrow, and their relationships are depicted in dark grey.
Relationships between taxa were determined following the phylogenies of the genera Vespa and Provespa (Saito & Kojima, 2011; Perrard et al., 2013) combined with the genetic analysis of V. velutina populations (Perrard et al., 2014; Fig. 2a). Two main genetic lineages were distinguished within V. velutina, but relationships within these lineages were not well supported and are considered here as unresolved. Species and population values used for the phylogenetic signal test were computed as the mean shape of each taxonomic group. Phy-
logenetic signal was tested with 10 000 permutations on both the raw and allometry-corrected data for the whole data set and the data set of Vespa specimens. Allometry, the influence of size on shapes, may have an influence on the phylogenetic signal because some species groups or genera have similar sizes, which may potentially induce a similarity of shapes within these groups. The evolutionary allometry (Klingenberg, 1996) has been estimated by a multivariate regression of the taxa mean wing shapes on the taxa mean wing sizes, defined as the
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
5
Results
whereas PCAPop, based on more specimens, accounted for only 6.03% (Table 2). The dimensionality of the variation in the tangent space increased progressively from the Vespinae genera subsample to the V. velutina populations’ subsample (Fig. 4). The first PC of PCAGen explained more of the total variation than the ones of PCASp and of PCAPop. In addition, fewer variables of PCAGen were required to describe 95% of the variation (Table 2). These results could not be explained by the smaller sample size and fewer of groups of PCAGen: the dimensionality increase remained clear even with subsamples of PCASp and of PCAPop having the same sample size and group number than PCAGen.
Variability of the wing shape in the tangent space
Gaussian mixture models
More than half of the variation in the Vespinae specimens was contained in one dimension. It described the variation between the rounded forewings of Provespa and the more elongated forewings of the other species (Figs 2b and 3). The wing shape variability decreased with the level of classification: variation in PCAGen reached 36.35% of the variation in the entire sample,
Gaussian mixture models applied on the entire data set resulted in 30 clusters highly congruent with species groups: 78.92% of the 1025 specimens were in the cluster of their own species, all specimens being in a cluster of their own genus. For 17 of the 28 species, every specimen of a same species was gathered in one cluster. For six species, including V. velutina, more than 15% of their specimens were included in a second cluster. These second groups were specific for two species. Vespa velutina specimens were grouped into six clusters, two of them being shared with Vespa vivax and Vespa bicolor specimens. In the subsample of genera, five Gaussian clusters corresponded to the four genera, Provespa specimens being grouped in two clusters. In the subsample of the 22 Vespa species, 18 clusters were found. Most of these clusters grouped individuals from one species, resulting in 85% of the 720 Vespa specimens grouped with individuals of their own species (Fig. 5). When applied to the 338 specimens of V. velutina, GMM resulted in only one cluster of all specimens.
log-transformed centroid size of the wing. Phylogenetic signal for corrected allometric data was computed on the residuals of this regression. Every analysis was conducted with the software R (R Development Core Team, 2014) using the packages ‘Rmorph’ (Baylac, 2013), ‘mclust’ (Fraley & Raftery, 2002), ‘MASS’ (Venables & Ripley, 2002) and ‘ape’ (Paradis et al., 2004) as well as original R functions for the test of the phylogenetic signal. These functions were based on the algorithm described by McArdle & Rodrigo (1994) for ancestral states reconstruction of continuous characters minimizing squared parsimony.
(a)
(b)
Differences between the taxonomic groups
(c)
Fig. 3 Wing shape variation between Vespa and Provespa. (a) Difference in mean shapes of Vespa (grey) and Provespa (black) wings. (b) Wing of Vespa luctuosa. (c) Wing of Provespa barthelemyi. Scale bar is 1 mm.
The CVA with cross-validation on the whole sample returned 88.29% of correct assignment of specimens to their species and populations. It confirmed that 100% of the 1025 specimens could be attributed to the right genus and 97.07% to the right species. Only 70.71% of V. velutina specimens were reassigned to the right population, but the result is quite unbalanced between continental and insular populations: 96.12% (125 of 129) of the insular specimens were attributed to their correct island, whereas only 55.02% (115 of 209) of the continental specimens were attributed to their own population. Correct assignment tended to increase with the number of variables and with the exclusion of small groups, but the trends between taxonomic levels remained similar (see Supporting Information). Reassignment based on the k-nearest neighbours also
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
6
A. PERRARD ET AL.
Table 2 Wing shape variation in the multivariate tangent space for the entire Vespinae sample and for each taxonomic level subsample.
Total sample Vespinae genera Vespa species Vespa velutina populations
N indiv.
N groups
r2 (10 5)
N PC (95%)
1st PC (%)
2nd PC (%)
3rd PC (%)
1025 159 720 338
36 4 22 9
126.0 45.8 40.8 7.6
12 6 16 18
54.68 69.94 30.35 27.26
12.59 13.46 19.07 14.05
7.02 6.36 11.46 9.49
40
60
80
PVespa = 0) and allometry-corrected data (D2total = 9.133 9 10 3, Ptotal = 0; D2Vespa = 4.071 9 10 3, PVespa = 0). The phylogeny mapped in the tangent space showed that wing shapes tended to be grouped according to their clade in most cases. In the subsample of V. velutina, the first axis distinguished populations of the two main genetic lineages. Nonetheless, the influence of the phylogeny on the variation in the wing venation was not complete, as shown by the many overlapping lineages and the diverging nodes in the tangent space. As such, strong divergences were found between the wing shapes of Provespa and those of the other vespine genera, between the Vespa mandarinia species group and the Vespa tropica species group, as well as between the three closely related species Vespa affinis, Vespa orientalis and Vespa mocsaryana (Fig. 2b,c).
Vespinae genera Vespa species V. velutina populations
20
Percentage of variation
100 95
The column ‘N indiv.’ indicates the sample size, ‘N groups’ refers to the number of taxonomic groups compared, ‘r2’ refers to the variance in the total shape space, and ‘N PC’ refers to the number of principal components (or independent dimensions) required to explain 95% of the shape variation in the tangent space comprising 34 dimensions. The last three columns describe the percentage of variation explained by the three-first axes of the principal component analyses for each sample. The variance in the total shape space was assessed by the sum of the PC scores’ variances normalized by the ratio between each subsample size and the total sample size.
0
5
10 15 20 25 Number of dimensions
30
34
Fig. 4 Cumulative percentage of the total variation explained by the first components of principal component analyses on the Vespinae genera, Vespa species and Vespa velutina populations’ subsamples. Error bars are 95% estimates from 1000 bootstrap replicates.
confirmed these CVA results. The variation in wing shape means was 6.43 P times higher between genera than P between species ( 2genera = 0.00274, Ngenera = 4 2 and species = 0.00045, Nspecies = 22). The variation between species was 4.54 times higher than P that found between populations of V. velutina ( 2populations = 0.00009, Npopulations = 9). However, this variation between populations of the same species was slightly lower than the overall intrapopulation variation estimated on population-centred shape data in V. velutina P ( 2intrapopulation = 0.00015, Nintrapopulation = 338). Phylogenetic signal Phylogenetic signal was significant both in the whole data set and in the Vespa subsample, with both raw (D2total = 9.144 9 10 3, Ptotal = 0; D2Vespa = 3.902 9 10 3,
Discussion Wing venation shape as a taxonomic marker Our results showed that the evolution of the wing venation may explain its discriminating power by inducing structured, original and progressive differences of the venation shape between the different taxa. The genetic drift between taxa may have thus played a major role in the wing evolution. First, the main structure of the wing shape variation was found congruent with the taxonomic groups of social wasps. Many components could have been involved in the variation of wing shapes of 1025 specimens originating from a wide range of localities and belonging to species with largely overlapping distributions (Carpenter & Kojima, 1997). Gaussian mixture models showed, nonetheless, that the distribution of wing shapes in the tangent shape-space reflected species or groups of closely related species. These results confirmed that, within sex and caste, the genetic difference between species had a stronger influence on the wing shape than other biotic or abiotic factors. Secondly, the wing shape variation was highly multivariate: its evolution did not follow simple variation patterns even if the variation was not high in all 34
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Fig. 5 Species composition of the 18 Gaussian clusters based on the wing shape distribution of the Vespa subsample. Most clusters group mainly specimens from one species or from closely related species. The largest cluster (N = 77) groups many specimens from various species without taxonomical tendency. Species legends are as follows: a. Vespa bicolor; b. Vespa crabro; c. Vespa fervida; d. Vespa luctuosa; e. Vespa mandarinia; f. Vespa mocsaryana; g. Vespa multimaculata; h. Vespa philippinensis; i. Vespa simillima; j. Vespa soror; k. Vespa vivax; l. Vespa analis; m. Vespa ducalis; n. Vespa bellicosa.
Vespa affinis V. analis V. basalis V. analis V. bellicosa V. bicolor V. simillima V. binghami V. crabro V. ducalis V. dybowskii V. fervida V. fumida V. mandarinia V. mocsaryana V. multimaculata V. orientalis V. tropica V. soror V. velutina 0.0
0.2
dimensions of the tangent space (Fig. 4). This high dimensionality of the shape variation was previously reported in other studies on flies and honey bees (Monteiro et al., 2002; Mezey & Houle, 2005). Under a model of random evolution, high dimensionality of variation increases the potential ways of divergence during the evolutionary history. As a result, once two groups have diverged, they are less likely to converge to an exact same location again in the multidimensional space (C.P. Klingenberg, personal communication). The relatively high dimensionality of the wing shape variation may have promoted the evolution of specific venation shapes by preventing complete convergence of the wing between different taxa. Furthermore, the wing shape variation not only reflected the taxonomic groups, but was also partially congruent with the evolutionary history shared by these groups. The CVA showed that the distribution of wing shapes of a single genus in the tangent space was not overlapping with the ones from its closest genera. On the other hand, the distinction was lower between species: 2.93% of the individuals were attributed to a wrong species within the same genus and even lower between populations of the single species V. velutina. The insular populations of V. velutina were more readily discriminated than the continental populations (96.12% vs. 55.02%). It is therefore possible that genetic isolation experienced in islands enhanced the evolution of distinct wing shapes (Camara et al., 2006). Furthermore, both the GMM and CVA misattributions between species mainly involved closely related species, suggesting that closely related species tended to have more similar wing shapes. The comparison of variances at the different levels of the classification confirmed that the variation between groups increased with the taxonomic level. Finally, the significance of the phylogenetic signal showed that the variation in wing shape
0.4
c. a. b.
7
n = 34 V. crabro n = 19 g. k. n = 27 V. velutina n = 77 f. h. i. j. V. velutina V. vivax n = 64 n = 24 n = 47 V. dybowskii n = 30 n = 20 n = 43 V. bellicosa n = 13 n = 22 n = 26 V. bicolor V. bellicosa l. n = 47 n = 38 V. philippinensis m. n. n = 43 n = 32 V. mandarinia n = 22
d. e.
V. luctuosa
0.6
0.8
1.0
between the populations, species and genera was congruent with their evolutionary history. The unresolved nodes in the phylogeny did not interfere with this result: by increasing the length of the tree in the tangent space (Klingenberg & Gidaszewski, 2010), unresolved nodes decrease the probability of a significant signal. The phylogenetic signal being based on a Brownian motion model of the phenotype evolution, this last result supports an evolution of the phenotype mainly driven by genetic drift. Altogether, these results are further evidence that the variation in the wing venation is largely shaped by the evolutionary history and could be used to identify the genera and species of insects presenting a venation pattern suitable for landmark analyses. Other factors of wing shape variation Our results also suggested that drift may not have been the only factor influencing the wing shape evolution in social wasps. Although there is a significant phylogenetic signal in the overall variation, some of our results were more likely explained by the influence of external factors than by pure drift.
Body shape The main wing shape variation found in our sampling was the difference between the apically rounded Provespa wings and the more elongated wings of other Vespinae (Figs 2b and 3). A similar variation of the wing was found in Polistinae where rounded wings were related to small body sizes (Garcıa & Sarmiento, 2012). In our case, this relationship does not hold as Provespa specimens were larger than those of Vespula and Dolichovespula. Some of the main differences between Provespa and other Vespinae genera such as enlarged ocelli and brown coloration are adaptations to nocturnal hab-
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
A. PERRARD ET AL.
(a)
0.01
1 v4 v1 2
v2 v3 v5
PC2 −0.01
0.00
3
v7 6
v6 4 v8 5
v9 7
−0.02
its. However, no straightforward link can be established between wing shape and nocturnal habits. Furthermore, one hornet species, Vespa binghami, is nocturnal and shares the enlarged ocelli and brown coloration with Provespa while retaining a Vespa-like wing shape. The divergent wings of Provespa species may thus be more likely related to their slender and proportionally lighter body than those of yellow jackets and hornets. The wing venation has a role in wing rigidity (Combes & Daniel, 2003) which is important for flight (Young et al., 2009; Mountcastle & Combes, 2013). Furthermore, the shape of the wing is likely involved in the production of the aerodynamic forces during the flight (Young et al., 2009). A variation in the ratio of the body mass compared with its length during the evolution of Provespa wasps may have induced a shift in the optimal wing shape. These divergent selective pressures could have enhanced the development of divergent wing shapes between Provespa species and the other Vespinae.
8
−0.04
8
−0.02
−0.01
0.01 0.00 PC1
0.02
0.03
(b) PC2
Climate The wing shape of V. orientalis, the only species of hornet living in arid climates, was distinguished from those of all other Vespa species by the second PC of variation within the genus (Fig. 6). Furthermore, the same axis of variation separated the mean shapes of continental and Indonesian populations of V. velutina, the latter being closer to the mean shapes of the four Indonesian species of the same clade than V. velutina. On the other hand, continental populations of V. velutina were closer to continental species of the clade including V. velutina (Fig. 6a). These results suggested that the second PC of the variation within Vespa described a shape change influenced by geographical and potentially climatic components. Latitudinal and altitudinal clines of wing variation were found in several species of Drosophila (e.g. Dahlgaard et al., 2001; Sambucetti et al., 2006; Pitchers et al., 2013). Furthermore, wing shape development is known to be influenced by temperature (Debat et al., 2006; Pitchers et al., 2013). The variation of climatic parameters such as humidity and temperature in different geographical ranges may have influenced the shape variation between wasp species either through plasticity of the phenotype or by producing different selective pressures. Phenotypic plasticity, the ability to develop different phenotypes with a same genetic background in response to different environments, has been previously found in wings of Drosophila in relation to thermal clines (Debat et al., 2006; Pitchers et al., 2013). Phenotypic plasticity is likely to have played a role in the observed variation of Vespa wings because this variation was detected between species as well as within a single species (Fig. 3). However, hornets develop in enclosed nests with buffered temperature which should limit the effect of climatic differences on wings develop-
Fig. 6 Variation in wing shape potentially influenced by climatic factors. (a) Two first dimensions of the tangent shape space of Vespa wings with the phylogeny mapped. The second principal component (PC) differentiates continental (white squares) and Indonesian (white circles) species and populations of the clade of Vespa velutina as well as the only species from arid zones Vespa orientalis (star). 1: Vespa vivax; 2: Vespa simillima; 3: Vespa bicolor; 4: Vespa multimaculata; 5: Vespa bellicosa; 6: Vespa luctuosa; 7: Vespa fervida; 8: V. orientalis; v1–9: V. velutina. v1: Vietnam; v2: Zhejiang; v3: Yunnan; v4: Thailand; v5: Nepal; v6: Java; v7: Sulawesi; v8: Flores; v9: Lombok. (b) Wing shape variation along the second PC axis from negative (grey) to positive (black) values.
ment (Jones & Oldroyd, 2006). Furthermore, the geographical clines could not be entirely explained by external factors even in Drosophila (Dahlgaard et al., 2001; Sambucetti et al., 2006). Part of the wing variation related to the geography in hornets, such as the highly divergent wings of V. orientalis, may thus come from fixed genetic differences driven by divergent selective pressures. Phenotypic plasticity promotes evolutionary answers to such selective pressures through ‘genetic accommodation’ (West-Eberhard, 2005). The combination of phenotypic plasticity and divergent selection may thus explain the observed variation of wing shape congruent with the geography.
Mimicry The wing shape of V. affinis, a species closely related to V. mocsaryana, had an extreme value on the first PC of the PCASp, close to the shape of its distant relative V. tropica (Fig. 2c). These two species present also similar geo-
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
graphical distribution and coloration, a dark body with a pale area on the two-first metasomal terga. On the other hand, wing shapes of V. mocsaryana and Vespa basalis were the closest, both species having reddish colour and almost superimposed distributional range (Carpenter & Kojima, 1997). The hypothesis of a convergence due to shared climatic conditions cannot be retained for these species as these similar distributions, and thus climatic areas, are largely overlapping with those of the other Vespa species. However, wing shape was related to wing sound production in Drosophila (Routtu et al., 2007; Menezes et al., 2013), and buzzing was proven to enhance the aposematism efficiency in wasps (Hauglund et al., 2006). Wing shape could thus be selected by mimicry pressures on the flight sound. A second explanation could be that genes influencing the patterning of pigmentation, such as genes Hox (Kopp et al., 1999; Perrard et al., 2014), may also be involved in the development of similar venation patterns. Further research is required to properly test the potential influence of these factors on wing shape and to understand the origin of the resulting variation, from plasticity or fixed genetic differences. Our results, nonetheless, demonstrate that these external influences are not strong enough to blur the wing shape differences due to the divergence between genera and species.
Conclusion Wing shape measurement and analysis is fast and requires only a reference data set, a camera and a computer. Previous studies have shown that this is an efficient identification criterion on a wide variety of insects. These include economically important species such as vectors of diseases (Jirakanjanakit & Dujardin, 2005; Henry et al., 2010), pests (Schutze et al., 2011) and pollinators (Francoy et al., 2008, 2009). Our results provide new insights confirming and explaining the efficiency of the venation shape for insect identification: the main source of variation in the venation shape was the divergence between taxa, and this divergence occurred progressively, with rare and incomplete secondary convergence. With an initial effort to produce the reference data sets, using identified specimens from natural history collections, for example, it could become a quick and inexpensive identification tool (Steinhage et al., 2001; Camara et al., 2006; Schutze et al., 2011). Furthermore, standardized wing imaging can be performed using portable material such as a camera with a tripod and microscope blades (Perrard et al., 2012). Wing shape could thus enable insect identification directly on the field, which would represent an important advance for many ecological studies and pest surveys. Finally, our results not only explained the efficiency of the wing shape as an identification criterion, but also showed clues to other factors influencing its evolution.
9
The wing venation has a clear functional role in flight by determining the rigidity of the wing (Combes & Daniel, 2003; Young et al., 2009; Mountcastle & Combes, 2013). However, the functional importance of the venation shape is still poorly understood (Routtu et al., 2007; Mountcastle & Combes, 2013). In our analysis, the wing shape variation tended to be highly multivariate and homogeneous within a species while shapes diverged along specific directions between distant species or genera. Whether these evolutionary pathways in the space of possible shapes resulted from developmental constraints or selective pressures is yet to be elucidated. Further analyses of these patterns of shape variation in regard to potential adaptive processes could help to understand the role of the wing shape and venation in insect flight.
Acknowledgments This work was funded by two ATM (‘Formes possibles, formes r ealis ees’ and ‘Biodiversit e actuelle et fossile’) grants from the Museum National d’Histoire Naturelle, the Bourse Eole of the Franco-Netherland Alliance and the Bourse Germaine Cousin of the Soci et e Entomologique de France. We thank C.P. Klingenberg, P. Mitteroecker and the anonymous reviewer for their helpful comments on the manuscript. We are grateful to J. Haxaire, J.-L. Renesson and P. Tripotin who shared their private collections. Special thanks are due to curators A. Drumont and G. Yvonnick (Institut Royal des Sciences Naturelles de Belgique, Brussels), G. Broad (Natural History Museum, London), F. and J. Gusenleitner (Ober€ osterreiches Landesmuseum, Linz), K. van Achterberg (Nationaal Natuurhistorisch Museum Naturalis, Leiden), H. Nugroho (Museum Zoologicum Bogoriense, Bogor) and X. Chen (University of Zhejiang, Hangzhou) who allowed to study in their respective museum collections. This study also benefited from the public collections of the American Museum of Natural History (New York, USA) and the Mus eum National d’Histoire Naturelle (MNHN, Paris, France). This study is a contribution of the Morphometric facility of the MNHN (plate-forme Morphom etrie, CNRS-MNHN UMS 2700).
Data accessibility Landmark coordinates of the measured wings, R scripts and R functions to compute the phylogenetic signals are available from the Dryad Digital Repository: http://doi. org/10.5061/dryad.2cv1h.
References Adams, D.C., Rohlf, F.J. & Slice, D.E. 2004. Geometric morphometrics: ten years of progress following the ‘revolution’. Ital. J. Zool. 71: 5–16.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
10
A. PERRARD ET AL.
Bai, M., McCullough, E., Song, K.-Q., Liu, W.-G. & Yang, X.-K. 2011. Evolutionary constraints in hind wing shape in Chinese dung beetles (Coleoptera: Scarabaeinae). PLoS ONE 6: e21600. Baylac, M. 2013. Rmorph: a R geometric and multivariate morphometrics library. Available from the author: [email protected] Baylac, M. & Daufresne, T. 1996. Wing venation variability in Monarthropalpus buxi (Diptera, Cecidomyiidae) and the quaternary coevolution of box (Buxus sempervirens L.) and its midge: a geometrical morphometric analysis. In: Advances in Morphometrics (L. Marcus, ed.), pp. 285–302. Plenum Press, New York, NY. Baylac, M., Villemant, C. & Simbolotti, G. 2003. Combining geometric morphometrics with pattern recognition for the investigation of species complexes. Biol. J. Linn. Soc. 80: 89– 98. Bickford, D., Lohman, D.J., Sodhi, N.S., Ng, P.K., Meier, R., Winker, K. et al. 2007. Cryptic species as a window on diversity and conservation. Trends Ecol. Evol. 22: 148–155. Bookstein, F.L. 1991. Morphometric Tools for Landmark Data: Geometry and Biology. Cambridge University Press, New York, NY. Camara, M., Caro-Ria~ no, H., Ravel, S., Dujardin, J., Hervouet, J., De Mee€ us, T. et al. 2006. Genetic and morphometric evidence for population isolation of Glossina palpalis gambiensis (Diptera: Glossinidae) on the Loos islands, Guinea. J. Med. Entomol. 43: 853–860. Carpenter, J.M. & Kojima, J. 1997. Checklist of the species in the subfamily Vespinae (Insecta: Hymenoptera: Vespidae). Nat. Hist. Bull. Ibaraki Univ. 1: 51–92. Combes, S.A. & Daniel, T.L. 2003. Flexural stiffness in insect wings I. Scaling and the influence of wing venation. J. Exp. Biol. 206: 2979–2987. Comstock, J.H. & Needham, J.G. 1899. The wings of insects. Am. Nat. 33: 845–860. Dahlgaard, J., Hasson, E. & Loeschcke, V. 2001. Behavioral differentiation in oviposition activity in Drosophila buzzatii from highland and lowland populations in Argentina: plasticity or thermal adaptation? Evolution 55: 738–747. De Meulemeester, T., Michez, D., Aytekin, A.M. & Danforth, B.N. 2012. Taxonomic affinity of halictid bee fossils (Hymenoptera: Anthophila) based on geometric morphometrics analyses of wing shape. J. Syst. Palaeontol. 10: 755–764. Debat, V., Milton, C.C., Rutherford, S., Klingenberg, C.P. & Hoffmann, A.A. 2006. Hsp90 and the quantitative variation of wing shape in Drosophila melanogaster. Evolution 60: 2529– 2538. Dryden, I.L. & Mardia, K.V. 1998. Statistical Shape Analysis. John Wiley & Sons, Chichester. Fraley, C. & Raftery, A.E. 1998. How many clusters? Which clustering method? Answers via model-based cluster analysis. Comput. J. 41: 578–588. Fraley, C. & Raftery, A.E. 2002. Model-based clustering, discriminant analysis, and density estimation. J. Am. Stat. Assoc. 97: 611–631. Francoy, T.M., Wittmann, D., Drauschke, M., M€ uller, S., Steinhage, V., Bezerra-Laure, M.A. et al. 2008. Identification of Africanized honey bees through wing morphometrics: two fast and efficient procedures. Apidologie 39: 488–494. Francoy, T.M., Silva, R.A.O., Nunes-Silva, P., Menezes, C. & Imperatriz-Fonseca, V.L. 2009. Gender identification of five
genera of stingless bees (Apidae, Meliponini) based on wing morphology. Genet. Mol. Res. 8: 207–214. Garcıa, Z. & Sarmiento, C.E. 2012. Relationship between body size and flying-related structures in neotropical social wasps (Polistinae, Vespidae, Hymenoptera). Zoomorphology 131: 1– 11. Gaston, K.J. & O’Neill, M.A. 2004. Automated species identification: why not? Philos. Trans. R. Soc. Lond. B Biol. Sci. 359: 655–667. Grimaldi, D.A. & Engel, M.S. 2005. Evolution of the Insects. Cambridge University Press, Cambridge. Hauglund, K., Hagen, S.B. & Lampe, H.M. 2006. Responses of domestic chicks (Gallus gallus domesticus) to multimodal aposematic signals. Behav. Ecol. 17: 392–398. Henry, A., Thongsripong, P., Fonseca-Gonzalez, I., Jaramillo-Ocampo, N. & Dujardin, J.-P. 2010. Wing shape of dengue vectors from around the world. Infect. Genet. Evol. 10: 207–214. Huey, R.B. & Bennett, A.F. 1987. Phylogenetic studies of coadaptation: preferred temperatures versus optimal performance temperatures of lizards. Evolution 41: 1098–1115. Jirakanjanakit, N. & Dujardin, J.-P. 2005. Discrimination of Aedes aegypti (Diptera: Culicidae) laboratory lines based on wing geometry. SE. Asian J. Trop. Med. 36: 858–861. Jones, J.C. & Oldroyd, B.P. 2006. Nest thermoregulation in social insects. Adv. In Insect Phys. 33: 153–191. Kendall, D.G. 1977. The diffusion of shape. Adv. Appl. Probab. 9: 428–430. Klingenberg, C.P. 1996. Multivariate allometry. In: Advances in Morphometrics (L. Marcus, ed.), pp. 23–49. Plenum Press, New York, NY. Klingenberg, C.P. & Gidaszewski, N.A. 2010. Testing and quantifying phylogenetic signals and homoplasy in morphometric data. Syst. Biol. 59: 245–261. Klingenberg, C.P., Barluenga, M. & Meyer, A. 2002. Shape analysis of symmetric structures: quantifying variation among individuals and asymmetry. Evolution 56: 1909–1920. Kopp, A., Blackman, R.K. & Duncan, I. 1999. Wingless, decapentaplegic and EGF receptor signaling pathways interact to specify dorso-ventral pattern in the adult abdomen of Drosophila. Development 126: 3495–3507. Marcus, L., Hingst-Zaher, E. & Zaher, H. 2000. Application of landmark morphometrics to skulls representing the orders of living mammals. Hystrix 11: 27–47. Matias, A., de la Riva, J., Torrez, M. & Dujardin, J.-P. 2001. Rhodnius robustus in Bolivia identified by its wings. Mem. I. Oswaldo Cruz 96: 947–950. Mayr, E. 1969. Principles of Systematic Zoology. McGraw-Hill, New York, NY. McArdle, B. & Rodrigo, A.G. 1994. Estimating the ancestral states of a continuous-valued character using squared-change parsimony: an analytical solution. Syst. Biol. 43: 573–578. Meier, R., Shiyang, K., Vaidya, G. & Ng, P.K. 2006. DNA barcoding and taxonomy in Diptera: a tale of high intraspecific variability and low identification success. Syst. Biol. 55: 715–728. Menezes, B.F., Vigoder, F.M., Peixoto, A.A., Varaldi, J. & Bitner-Mathe, B.C. 2013. The influence of male wing shape on mating success in Drosophila melanogaster. Anim. Behav. 85: 1217–1223. Mezey, J.G. & Houle, D. 2005. The dimensionality of genetic variation for wing shape in Drosophila melanogaster. Evolution 59: 1027–1038.
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
Wing shape evolution and species recognition
Mitteroecker, P. & Bookstein, F. 2011. Linear discrimination, ordination, and the visualization of selection gradients in modern morphometrics. Evol. Biol. 38: 100–114. Monteiro, L.R., Diniz-Filho, J.A.F., Reis, S.F. & Ara ujo, E.D. 2002. Geometric estimates of heritability in biological shape. Evolution 56: 563–572. Mountcastle, A.M. & Combes, S.A. 2013. Wing flexibility enhances load-lifting capacity in bumblebees. Proc. Roy. Soc. B: Biol. Sci. 280: 20130531. Outomuro, D., Adams, D.C. & Johansson, F. 2013. The evolution of wing shape in ornamented-winged damselflies (Calopterygidae, Odonata). Evol. Biol. 40: 300–309. Paradis, E., Claude, J. & Strimmer, K. 2004. APE: analyses of phylogenetics and evolution in R language. Bioinformatics 20: 289–290. Perrard, A., Villemant, C., Carpenter, J.M. & Baylac, M. 2012. Differences in caste dimorphism among three hornet species (Hymenoptera: Vespidae): forewing size, shape and allometry. J. Evol. Biol. 25: 1389–1398. Perrard, A., Pickett, K.M., Villemant, C., Kojima, J. & Carpenter, J.M. 2013. Phylogeny of hornets: a total evidence approach (Hymenoptera, Vespidae, Vespinae, Vespa). J. Hymenopt. Res. 32: 1–15. Perrard, A., Arca, M., Rome, Q., Muller, F., Tan, J., Bista, S. et al. 2014. Geographic variation of melanisation patterns in a hornet species: genetic differences, climatic pressures or aposematic constraints? PLoS ONE 9: e94162. Pitchers, W., Pool, J.E. & Dworkin, I. 2013. Altitudinal clinal variation in wing size and shape in African Drosophila melanogaster: one cline or many? Evolution 67: 438–452. R Development Core Team. 2014. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. http://R-project.org (accessed 16 March 2014). Ripley, B.D. 1996. Pattern Recognition and Neural Networks. Cambridge University Press, Cambridge. Rohlf, F.J. 2010. tpsDig v2. 16. Department of Ecology and Evolution, State University of New York, Stony Brook, NY (Distributed by author). Routtu, J., Mazzi, D., Van Der Linde, K., Mirol, P., Butlin, R.K. & Hoikkala, A. 2007. The extent of variation in male song, wing and genital characters among allopatric Drosophila montana populations. J. Evol. Biol. 20: 1591–1601. Saito, F. & Kojima, J. 2011. Phylogenetic analysis and biogeography of the nocturnal hornets, Provespa (Insecta: Hymenoptera: Vespidae: Vespinae). Species Div. 16: 65–74. Sambucetti, P., Loeschcke, V. & Norry, F.M. 2006. Developmental time and size-related traits in Drosophila buzzatii along an altitudinal gradient from Argentina. Hereditas 143: 77–83.
11
Schutze, M.K., Jessup, A. & Clarke, A.R. 2011. Wing shape as a potential discriminator of morphologically similar pest taxa within the Bactrocera dorsalis species complex (Diptera: Tephritidae). Bull. Entomol. Res. 102: 103–111. Steinhage, V., Arbuckle, T., Schr€ oder, S., Cremers, A.B. & Wittmann, D. 2001. ABIS: automated identification of bee species. In: BIOLOG. Workshop, German Programme on Biodiversity and Global Change Status Report (W.G. Berendsohn, ed.), pp. 194–195. Bonn. Venables, W.N. & Ripley, B.D. 2002. Modern Applied Statistics with S, 4th edn. Springer Verlag, Berlin. Villemant, C., Simbolotti, G. & Kenis, M. 2007. Discrimination of Eubazus (Hymenoptera, Braconidae) sibling species using geometric morphometrics analysis of wing venation. Syst. Entomol. 32: 625–634. West-Eberhard, M.J. 2005. Developmental plasticity and the origin of species differences. Proc. Natl. Acad. Sci. USA 102: 6543–6549. Young, J., Walker, S.M., Bomphrey, R.J., Taylor, G.K. & Thomas, A.L.R. 2009. Details of insect wing design and deformation enhance aerodynamic function and flight efficiency. Science 325: 1549–1552. Zelditch, M., Swiderski, D.L., Sheets, H.D. & Fink, W.L. 2004. Geometric Morphometrics for Biologists: A Primer. Elsevier Academic Press, London.
Supporting information Additional Supporting Information may be found in the online version of this article: Table S1 Classification to genera. Table S2 Classification to species. Table S3 Classification to species excluding groups of sample sizes < 10. Table S4 Classification to species excluding groups of sample sizes < 15. Table S5 Classification to species excluding groups of sample sizes < 20. Table S6 Classification to populations. Table S7 Classification to populations excluding groups of sample size < 10. Table S8 Classification to populations excluding groups of sample size < 34. Received 29 May 2014; revised 18 September 2014; accepted 26 September 2014
ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY. J. EVOL. BIOL. doi: 10.1111/jeb.12523 JOURNAL OF EVOLUTIONARY BIOLOGY ª 2014 EUROPEAN SOCIETY FOR EVOLUTIONARY BIOLOGY
