Butterflies and moths exhibit a spectacular diver sity of wing shape and size. The extent of wing variation is particularly evident in wild silk moths (Saturniidae), which have large wing shape and size variation. Some species have jagged wing margins, rounded forewing apical lobes, or narrow hind wings with long tails, while others lack these traits entirely. Surprisingly, very little work has been done to formally quantify wing variation within the fa mily. We analyzed the hind wing shape and size of 76 saturniid species representing 52 genera across five subfamilies using geometric morphometrics. We identified fifteen landmarks that we predict can be applied to families across Lepidoptera. PCA analyses grouped saturniid hind wings into six distinct morphological clusters. These groups did not appear to follow species relatedness—some phylogenetically and genetically distantly related taxa clustered in the same morphological group. We discuss ecological factors that might have led to the extraordinary wing variation within Saturniidae.
Butterflies and moths are known for their incredible diversity in wing morphology, which has fascinated naturalists for centuries (e.g., Merian 1705, Bates 1862, Grote 1897, Field 1898, Howse & Wolfe 2012). Some families have wings that resemble plumes (Alucitidae, Pterophoridae), while others are wasp-like (many Sesiidae). Perhaps the greatest degree of wing variation within a single lepidopteran family occurs in the Saturniidae, a cosmopolitan moth family that includes over 2,300 described species (van Nieukerken et al. 2011). Saturniidae contains some of the world's largest insects, with species that have wingspans of up to 250 mm, while others are small, with wingspans of only 30 mm (Imes 1992). Saturniids vary greatly in wing shape; some species have leaf-like jagged wing margins, expanded round forewing apical lobes, or have narrow hind wings with long tails, while others lack them entirely (Ylla et al. 2005, Nath & Devi 2009). Thus far, research on wing variation has focused primarily on venation (e.g., Packard 1905, Michener 1952, Albrecht & Kaila 1997) or on color pattern (e.g., Beldade & Brakefield 2002, Prieto et al. 2009, Monteiro et al. 1994, Mallet & Gilbert 1995). The few studies on moth wing shape have primarily dealt with intraspecific variation and sexual dimorphism (Ricklefs & O'Rourke 1975, Ricklefs 2009, Benitez et al. 2011, Nath & Devi 2009).
In Lepidoptera, wing morphology can be correlated with dispersal (Hill et al. 1999, Altizer & Davis 2010, Sekar 2012), gliding (Cespedes et al. 2015), mating (Wickman et al. 1992), predator evasion (Barber et al. 2015), or possibly nectar searching (Betts & Wootton 1988). However, studies that correlate wing size and wing shape with a species' biology across many taxa are limited, in part because comprehensive wing morphology data and robust phylogenies are often lacking. In this study, we characterized hind wing morphology in 76 distantly related saturniid species using geometric morphometrics and inferred the evolution of wing morphology using genetic similarity and a phylogeny. We focused on the hind wing, as it is an area of the wing that is highly variable and has direct implications for ecology, especially predator-prey interactionsBarber et al. 2015).
Materials and Methods
Specimen selection
We sampled a phylogenetically and morphologically diverse set of 76 saturniid species, representing 52 genera in five subfamilies across a broad geographical range (Supplementary Table 1). Specimens were obtained from the dry, pinned collection of the McGuire Center for Lepidoptera and Biodiversity, Florida Museum of Natural History (MGCL), Gainesville, Florida, USA. All specimens studied had undamaged right hind wings. In order to avoid bias associated with sexual dimorphism, we examined males. Most species were represented by a single specimen; five species were represented by ≥ 6 specimens as exemplars to assess the degree of intraspecific variation (Supplementary Table 1). Statistically significant low intraspecific hind wing variation (see Results) justified the examination of one specimen per species. All specimens chosen for inclusion in the study were digitally imaged for geometric morphometrics, and these pinned specimens have an associated white label with the printed text, “tail project.” Supplementary tables and additional files can be found on the Dryad Data Repository ( www.datadryad.org, accession number 10.5061/dryad.gs296).
Digitization
In order to obtain high-quality images for digitization, each saturniid specimen was pinned to an upright, square, 355 × 355 × 3 mm opaque acrylic sheet placed 18 cm above a millimeter ruler, positioned parallel to the focal plane of the right hind wing. The pin on the moth was inserted into a small foam block glued to one side of the sheet. A digital image of each specimen was taken in RAW format (5184 × 3456 pixels) using an EOS Rebel T3i digital camera with EF 35 mm F2.0 lens, mounted on an aluminum tripod 50 cm from the acrylic sheet. Two external Yongnuo YN 560 III flashguns, one placed behind the acrylic sheet, and another in front of it, were positioned to overexpose wings so that the venation could be visualized clearly. Prior to imaging, the dorsal surface of the right hind wing was lightly sprayed with 70% ethanol. Digital images of each specimen were taken in RAW format, converted to TIFF, and labeled in Photoshop CS5® with a species name and a unique identification number (e.g., “LEP01234“).
Tps file compilation and landmark designation in tpsDig
Each TIFF file was imported into tpsUtil ver. 1.58 (Schutz & Krieger 2007) and fifteen points on the dorsal surface of each hind wing were selected as homologous landmarks. Each wing image in the tps file was landmarked with tpsDig 2.12 (Rohlf 2010). Landmarks were chosen based on the approach of Zelditch et al. (2012); they were selected given the following criteria: a point was considered a good landmark if it was 1) easy to identify, 2) homologous and independent from other landmarks, 3) situated on the same plane as other landmarks, and 4) described the overall morphological shape of the wing being measured. Landmarks were located either at a point where a vein met the wing margin or where veins intersected each other (Fig. 1).
In order to accurately characterize wing shape, semilandmarks were also selected to supplement homologous landmarks. One hundred semilandmarks were placed along the perimeter of the right hind wing using tpsDig. The semilandmarks were applied in a clockwise direction from the base of the wing. Because tpsDig requires a consistent number of semilandmarks across species, the number of semilandmarks chosen was determined by the taxon that required the largest number of semilandmarks (i.e., Copiopteryx semiramis).
Procrustes Analysis, Principal Component Analysis, and estimating wing size
We used a Procrustes Analysis (PA) on the fifteen homologous landmarks to determine how each landmark varied among species and to determine if landmarks varied in a correlated fashion. In PA, landmarks are superimposed to optimize the squared distances from a common centroid through translation, rotation, and scaling, resulting in new Procrustes shape coordinates (Gower 1975). These coordinates are reduced to a smaller number of variables using a standard Principal Component Analysis (PCA). The Procrustes fit operation corrects for landmark errors by separating size and shape, and also allows shape coordinates to be projected into a Euclidean space tangent to the Procrustes shape space (Viscosi & Cardini 2012). After performing the PA, we retained only the axes with eigenvalues greater than 1 to represent the fifteen original landmarks. These analyses were performed in MorphoJ (Klingenberg 2011) and R (R Core Team 2014), using the R package, Shapes (Dryden 2012). The R script used for this analysis is provided in Appendix 1.
Fig. 1.
Digital images of the dorsal surface of the right hind wing of (A) Argema mimosae (Boisduval) (Saturniinae), and (B) Caio championi (Druce) (Arsenurinae). Numbers indicate the 15 homologous points chosen as landmarks. Landmark 1 was placed at the wing base, landmarks 2–9 were placed where Sc+R1 Rs M1 M2, M3, Cu1 Cu2, and A1+2 meet the wing margin, and landmarks 10–15 were placed around the periphery of the discal cell.
Cluster analysis
Using the results of the PA, we performed a cluster analysis to determine the number of groups that reduced the sum of squares of the data. To determine the most likely number of clusters in which the data could be split, we used the gap statistic as a measurement of goodness of clustering using fifteen as the maximum number of clusters (Tibshirani et al. 2001). Statistical significance was determined by the overlap between the confidence intervals estimated from 1000 bootstrap replicates. The cluster analysis was performed using the R package, Cluster (Maechler et al. 2013).
Estimating wing shape and size with semilandmarks
We performed an additional PA using the 100 semilandmarks in order to compare the results to those from the landmark data PA. Variation in wing shape was assessed using Procrustes ANOVA, using the groups from the landmark data PA as an independent variable. The Procrustes ANOVA calculates the sum of squares of the Procrustes distances between individuals, determines whether this sum is reduced by the introduction of an explanatory variable, and assesses statistical significance by randomizing individuals per group 999 times (Anderson 2001). Centroid size per individual, as estimated by the PA, was used as a proxy for wing size. Centroid size data was log-transformed to meet the assumptions for normality and one-way ANOVA was used to evaluate differences in wing size. All of the semilandmark analyses were performed using the R package, Geomorph (Adams & Otarola-Castillo 2013).
Estimating the evolution of wing morphology: phylogeny and genetic distance
In order to infer the evolution of hind wing morphology, we took two approaches: 1) phylogeny and 2) genetic distance. We used two methods because a comprehensive phylogeny including all 76 species examined in this project did not exist at the time of this study. To acquire a general understanding of how the six morphological groups evolved, we mapped the six groups onto a recently published phylogeny of Saturniidae (Barber et al. 2015). However, because only 12 of the 76 species studied were represented in the Barber et al. tree, we reduced it to a genera-level phylogeny and examined broad patterns of wing morphological evolution. This approach assumes that morphological groups do not vary within genera.
We also applied a genetic distance-based approach because it is conceivable that a genus includes species with different morphological groups. In order to determine if genetic distance correlates with morphological distance, we used the 658 bp, first subunit of the mitochondrial cytochrome oxidase gene (COI). COI is frequently used to delimit species of Lepidoptera, as there is a large reference sequence database of species (Hajibabaei et al. 2006, Hebert et al. 2003, Janzen et al. 2005). We downloaded all publicly available COI sequences from GenBank and the BOLD Taxonomy Browser (Ratnasingham & Hebert 2007) for the taxa in our study (Supplementary Table 1). We aligned the sequences using ClustalW (Larkin et al. 2007), implemented in Geneious v. 7.1 (Biomatters, Auckland, New Zealand) and estimated the best model of evolution using jModelTest v. 2.1.7 (Darriba et al. 2012). W used the R package, ape (Paradis et al. 2004) to estimate the corrected distance among species and performed a Principal Coordinates Analysis (PCoA) to scale the distance matrix into a twoaxis ordination. We used the PCoA scores to perform an analysis of similarity (ANOSIM) on the morphological groups, using the R package, Vegan (Oksanen et al. 2013). The ANOSIM enabled us to determine whether the calculated genetic distance between morphological groups is greater than expected by chance.
Finally, because we did not have a complete sampling for molecular data across all species, we calculated the genus to species ratio (Generic coefficient, Jaccard 1922) in each group. Traditionally, this ratio was used originally to represent the ecological diversity in a community, assuming that higher the ratio, the more ecologically diverse the community was (Jaccard 1922, Elton 1946). Even though some statistical sampling issues have been raised against the use of the genus to species ratio (Jarvinen 1982), if coupled with the appropriate null models, these ratios can be useful for the purposes of determining the diversity of a sample (Jarvinen 1982, Gotelli 2001). In our case, we did not use the genus to species ratio to infer the ecological diversity of a community but rather applied the same approach to represent phylogenetic diversity of morphological groups. Assuming that taxonomy is a fair representation of phylogenetic relationships, this method can be useful to make inferences about trait evolution in taxa with little or no molecular data.
In this approach, the genus to species ratio can vary from 1 to 1/n, where n is the group size. A value of 1 indicates that all species present in the group are representatives of different genera and it takes the value of 1/n when all species in the group belong to a single genus. To estimate the significance of the genus to species ratio, we randomly assigned a morphological group (1–6) to each species, constraining the number of species in each group to match the observed number of species per group in our data set. We then calculated the genus to species ratio for the randomly generated data set. We repeated this procedure 9999 times to create a null distribution of genus to species ratios to which we compared our observed values (Gotelli 2001). Using the null distribution, we calculated the standardized effect size of the genus to species ratio:
where G/S obs is the observed generic coefficient, G/Snull is the mean generic coefficient from the randomly generated groups, and G/S schnull is the standard deviation of the null distribution of generic coefficients. Values greater than 1 indicate that the generic coefficient is larger than expected by chance. Values smaller than 1 indicate that the generic coefficient is smaller than expected by chance. The P value for the genus to species ratio and standardized effect size was calculated using the formula:
Intraspecific versus interspecific variation
Since most of our sampling included only one individual per species, it was important to determine whether sampling limitation introduced a bias that would skew our results. We used a one-way ANOVA to compare the log-transformed variance of wing shape and size variables between individuals in the same species with the log-transformed variance of individuals belonging to different genera. However, because some morphological groups contained species belonging to different genera, we also compared the intraspecific variation with the intra-group variation. The objective of these comparisons is to demonstrate that the individuals sampled represent the mean value of wing shape and size for that species, and that intraspecific variation does not influence the patterns observed (Harmon & Losos 2005).
Fig. 2.
PCA results from Morphoj representing (A) PC1 and (B) PC2. The length and direction of the line corresponding to each point explains the variation found within each landmark and whether it is negatively or positively correlated to that PC. The x and y axes represent arbitrary coordinate values from a centroid in the original specimen image.
Results
Procrustes and Principal Components Analysis
The PA and ensuing PCA produced two PCs with eigenvalues greater than one, which explained 80% of the original variation. The first PC (PC1) was negatively correlated to the y values of landmarks 5 and 6 and positively correlated to the x and y values of 9 (Fig. 2A). The second PC (PC2) was negatively correlated to the y values of landmarks 2, 3, and 4, and positively correlated to the y values of 7 and the x and y values of 9 (Fig. 2B). All other landmarks were related to PCs with eigenvalues less than 1. The cluster analysis suggested that the species evaluated can be separated into six clusters (Supplementary Table 1).
Overall, saturniid hind wing morphology varied from relatively small, rounded wings (Group 1) to large wings with long, thin tails (Group 6, Supplementary Table 1). Species in Groups 1 and 2 had no hind wing tails. Groups 3 and 4 included species that have short, wide hind wings that either have very short tails or nubs. Group 3 (Arsenura, Caio) had hind wings that are more rounded than Group 4 (Dysdaemonia, Paradaemonia, Titaea). Group 5 was comprised of Actias and Graellsia that have shorter hind wing tails than those in Group 6. Group 6 consisted entirely of species in Copiopteryx.
Fig. 3.
Variation in hind wing shape across Saturniidae. (A) PCA results showing the six morphological groups based on the 15 homologous landmarks. Groups defined in this study are color-coded and the numbers at the top and right margins represent landmarks that are most highly correlated with PC1 and PC2, respectively. Black squares represent the centroid of each morphological group. Semi-transparent circles represent multiple individuals of the same species. (B) A comparison of COI genetic distance with the six morphological groups using multidimensional scaling. (C) A reduced, genus-level phylogeny, redrawn from Figure 2 of Barber et al. (2015) showing the six morphological groups.
PC1 characterized wing morphology based on tail length, whereas PC2 discriminated mainly on the width of the hind wing. PC1 thus created a gradient in which the tailless or very short tailed groups (Groups 1–4) were relatively close to each other (Fig. 3A). Although Groups 5 and 6 included tailed taxa, their PC groups were relatively distant from each other (Fig. 3A), suggesting that a tail is formed with either 1) having high x values of landmark 9 and intermediate y values of landmarks 7 and 9 (i.e. Group 5), or 2) having high x values of landmark 9 and high y values of landmarks 2, 3, 4, and 9 (i.e. Group 6).
Wing size and shape comparisons
The saturniids examined differed significantly in shape (Procrustes ANOVA: SS.obs = 2.97, df = 5, p = 0.001), with shapes following the pattern described using the PCA on the landmarks. Shape variation occurred on a gradient from hind wings that are rounded and tailless (Group 1) to hind wings that are elongate and tailed (Group 6; Fig. 4A). The species sampled also exhibited significant differences in wing size (ANOVA landmarks: F = 93.89, df = 5, p <0.001; semilandmarks: F = 13.68, df = 5, p<0.001). Size variation occurs on a gradient from small (Group 1) to large (Group 6; Fig. 4B). Interestingly, outlining the wing shape reveals that the hind wing tail is located in a similar position on most saturniids, but different wing veins form tails in unrelated groups.
Genetic distance and phylogeny
The COI genetic distance analysis revealed that the genetic distance within morphological groups is larger than the distance between the six morphological groups (ANOSIM: R statistic = -0.16, p = 0.049; Fig. 3B). Tailless groups were mainly composed of monotypic species, whereas tailed groups contain species that have multiple congeners in the same group (i.e., multiple species per genus; genus (G) to species (S) ratio: Group 1: G/S obs = 1.0, G/Sses = -1.38, p = 0.17; Group 2: G/Sobs = 1.0, G/Sses = -3.58, p < 0.01; Group 3: G/Sobs = 0.25, G/Sses= 7.1, p <0.01; Group 4: G/Sobs = 0.25, G/Sses = 7.05% < 0.01; Group 5: G/Sobs = 0.58, G/Sses = 3.81,p < 0.01; Group 6: G/Sobs = 0.17, G/Sses = 8.1 p < 0.01). Mapping the six groups onto the reduced, genus-level tree revealed that two morphological wing groups (Group 1 and 2) are present in multiple saturniid clades (Fig. 3C).
Intraspecific versus interspecific variation
Within each group, there was less intraspecific variation than interspecific variation (ANOVA; PC1: F = 10.75, df = 1, p < 0.01; PC2: t = 9.61, df = 1, p = 0.01). Hind wing variation within a species was equal to the variation within a single genus in PC1, but greater than the variation within a genus in PC2 (ANOVA; PC1: F = 3.7, df = 1, p = 0.08; PC2: t = 5.8, df = 1, p = 0.03), supporting the examination of one individual per species.
Discussion
Saturniid moths have evolved a broad range of hind wing shape and size. We identified fifteen points on the moth hind wing that can be consistently used as landmarks for geometric morphometrics. Of these landmarks, those positioned at the wing margin (landmarks 2–9; Fig. 2) exhibited the greatest interspecific variation. Landmarks located at margins of eyespots and wing marks (e.g., landmarks 12, 13; Fig. 1) were useful for intraspecific comparisons, but our tests for intraspecific variation were limited in sample size, and therefore should not be treated as final or conclusive. While saturniid hind wings were clustered into six distinct morphological groups, these groups did not correspond directly to COI genetic clusters. Mapping morphological groups onto phylogeny showed that wing shape and size might have evolved convergently among subfamilies of Saturniidae.
A question that remains largely unanswered is why unrelated saturniid species of different genera evolved significantly different wing shape and size. Wings of Lepidoptera are complex, serving multiple adaptive functions, such as flight, thermal regulation, mating behavior, and predator avoidance (reviewed in Scoble 1992). Hind wings have been shown to play a limited role in lepidopteran flight (Jantzen & Eisner 2008), but they might be critical under particular circumstances, such as evading deterring, or delaying predation (Collins 2013, Vallin et al. 2010).
Hind wing tails, which appear to have originated independently in several saturniid subfamilies, is thought to be an anti-bat strategy. Janzen (1984) postulated that moths with long tails make them appear larger to echolocating bats. Barber et al. (2015) showed that twisted hind wing tails are acoustic lures that reflect bat sonar and divert bat attack to these appendages. Insectivorous bats present a strong selective force on nocturnal moths (Conner & Corcoran 2012), and if the presence of a tail is an antibat defense, one could predict an evolutionary transition from an ancestral condition of no tails to long tails. This trend is seen across Saturniidae (Barber et al. 2015), and also within particular saturniid clades (e.g., the Actias group, Ylla et al. 2005), although these transitions are likely influenced by adaptive factors such as aerodynamic drag that could slow or change a moth's movement. A transition from simple anti-bat strategies to more complex ones is seen in other moth groups, such as hawkmoths, in which derived lineages use a combination of ultrasound hearing and jamming (Kawahara & Barber 2015). While we are beginning to understand behavioral and morphological adaptations of moth wings in light of evolution, what we know about them is still severely limited. Our next step is to examine moth wing traits (wing shape, size, color, pattern) and relate them to ecology and gene function.
Acknowledgements
We thank Chris Johns and Lary Reeves for providing digitization equipment and Ivonne Garzon for supplying morphometrics resources. Dean Adams assisted with the Geomorph package. James Adams, Jesse Breinholt, Jon Bremer, Jon Colburn, Michael Collins, Jaret Daniels, Dale Halbritter, Nick Homziak, Peter Houlihan, Simon McClung, Elena Ortiz-Acevedo, Pablo Sebastian Padron, Craig Segebarth, Ian Segebarth, John Snyder, Andrei Sourakov, Matt Standridge, Lisa Taylor, Andrew Warren and Lei Xiao provided support and feedback. This project was funded by the Mu Alpha Theta Summer high school research award to MZ, and in part by National Science Foundation grants DEB-15057007 and IOS-1121739 to AYK, IOS-1121807 to JRB, REU, ROA, and RAHSS supplements to these NSF awards.
