Evolution and losses of spines in slug caterpillars (Lepidoptera: Limacodidae)

Abstract Larvae of the cosmopolitan family Limacodidae, commonly known as “slug” caterpillars, are well known because of the widespread occurrence of spines with urticating properties, a morpho‐chemical adaptive trait that has been demonstrated to protect the larvae from natural enemies. However, while most species are armed with rows of spines (“nettle” caterpillars), slug caterpillars are morphologically diverse with some species lacking spines and thus are nonstinging. It has been demonstrated that the evolution of spines in slug caterpillars may have a single origin and that this trait is possibly derived from nonstinging slug caterpillars, but these conclusions were based on limited sampling of mainly New World taxa; thus, the evolution of spines and other traits within the family remains unresolved. Here, we analyze morphological variation in slug caterpillars within an evolutionary framework to determine character evolution of spines with samples from Asia, Australia, North America, and South America. The phylogeny of the Limacodidae was reconstructed based on a multigene dataset comprising five molecular markers (5.6 Kbp: COI, 28S, 18S, EF‐1α, and wingless) representing 45 species from 40 genera and eight outgroups. Based on this phylogeny, we infer that limacodids evolved from a common ancestor in which the larval type possessed spines, and then slug caterpillars without spines evolved independently multiple times in different continents. While larvae with spines are well adapted to avoiding generalist predators, our results imply that larvae without spines may be suited to different ecological niches. Systematic relationships of our dataset indicate six major lineages, several of which have not previously been identified.

the evolution of spines in slug caterpillars may have a single origin and that this trait is possibly derived from nonstinging slug caterpillars, but these conclusions were based on limited sampling of mainly New World taxa; thus, the evolution of spines and other traits within the family remains unresolved. Here, we analyze morphological variation in slug caterpillars within an evolutionary framework to determine character evolution of spines with samples from Asia, Australia, North America, and South America.
The phylogeny of the Limacodidae was reconstructed based on a multigene dataset comprising five molecular markers (5.6 Kbp: COI, 28S, 18S, EF-1α, and wingless) representing 45 species from 40 genera and eight outgroups. Based on this phylogeny, we infer that limacodids evolved from a common ancestor in which the larval type possessed spines, and then slug caterpillars without spines evolved independently multiple times in different continents. While larvae with spines are well adapted to avoiding generalist predators, our results imply that larvae without spines may be suited to different ecological niches. Systematic relationships of our dataset indicate six major lineages, several of which have not previously been identified.
Antipredator strategies occur in every biome of the world, implying that predation is a potent selective force and thus of immense ecological and evolutionary significance (Grimaldi & Engel, 2005;Murphy, Leahy, Williams, & Lill, 2010;Ruxton, Sherratt, & Speed, 2004). Spines are one kind of obvious antipredator strategy to avoid predation (Inbar & Lev-Yadun, 2005), such as the spines on inflated pufferfish (Brainerd, 1994), sticklebacks (Gross, 1978;Hoogland, Morris, & Tinbergen, 1956;Reimchen, 1983), slug caterpillars of the moth family Limacodidae (Murphy et al., 2010) and those on spiny plants (Gowda, 1996;Hanley, Lamont, Fairbanks, & Rafferty, 2007;Lev-Yadun, 2001). Spines are a common defense mechanism that have evolved independently (homoplasy) in aquatic and terrestrial ecosystems, indicating that the reappearance of this phenotype is highly adaptive. However, antipredator strategies may be secondarily lost due to various factors, for example, due to the loss of predators or limited nutrients (Bell, Francis, & Havens, 1985;Giles, 1983;Larson, 1976;McNab, 1994;Whitwell et al., 2012). Thus, it may be difficult to distinguish whether similar phenotypes present in a broadly distributed taxonomic clade is due to gains or losses. Hence, integrating phenotypic variation and reconstructing the probable ancestral states within a phylogenetic framework can enhance our knowledge of how traits evolve and may provide insights into the evolutionary processes and selective pressures involved.
Nettle caterpillars and gelatine caterpillars are almost distributed globally, whereas monkey slugs are rare, occurring in low abundance and being geographically restricted to Asia and the New F I G U R E 1 Different larval types of slug caterpillars in the Limacodidae with respect to the presence of spines: (a-c) first, early, and late instar of Parasa consocia (character state A: spines present after second instar); (d) late instar of Microleon longipalpis (character state A); (e-f) first and late instar larva of Cania heppneri (character state A); (g) spines on the late instar of Cania heppneri; (h) spines on the late instar of Microleon longipalpis; (i) first instar of Demonarosa rufotessellata subrosea (character state B: spines present after second instar but reduced in late instars); (j) second instar of Demonarosa rufotessellata subrosea with spines on the segments (character state B); (k) late instar of Demonarosa rufotessellata subrosea with almost all spines lost (character state B); (l) first instar of Phrixolepia inouei (character state D: spines absent but numerous setae present after second instar); (m) first instar of Caiella pygmy (character state B); (n) early instar of Caiella pygmy with spines (character state B); (o) late instar of Caiella pygmy with almost all spines reduced (character state B); (p) late instar of Phrixolepia inouei with numerous setae (character state D); (q) first instar of Pseudanapaea transvestita (character state B); (r) second instar of Pseudanapaea transvestita with spines (character state B); (s) late instar of Pseudanapaea transvestita with almost all spines reduced (character state B); (t) late instar larva of Nagodopsis shirakiana (character state C: spines absent in all instars); (u) early instar of Ecnomoctena brachyopa with spines (character state B); (v) late instar of Ecnomoctena brachyopa with almost all spines reduced (character state B); (w, x) first and late instar of Altha melanopsis (character state C) World. The majority of limacodid larvae are nettle caterpillars, which are armed with spines that are well known to inflict stings (Hossler, 2010;Kawamoto, 1978;Murphy et al., 2010;Walker, 2018;Zaspel et al., 2016). Murphy et al. (2010) presented evidence that spines do indeed protect slug caterpillars from generalist predators. Cock et al. (1987) presented a hypothesis that nonstinging types of slug caterpillars evolved from nettle caterpillars. However, the first detailed phylogenetic study of Limacodidae by Zaspel et al. (2016) suggested that (a) nettle caterpillars are a monophyletic group; (b) gelatine caterpillars are a monophyletic group; and (c) nettle caterpillars are derived from gelatine caterpillars. Because the study of Zaspel et al. (2016) was based on mainly New World taxa, the results may be derived from in situ diversification or independent colonization. Thus, it is uncertain if the evolutionary pattern of slug caterpillars is the same after including samples from different zoogeographic regions of the world. It is also unclear whether the existence or loss of spines in slug caterpillars has evolved once or has evolved repeatedly and independently in different lineages and/or in different continents.
When spines are present, they may be derived from a common ancestor or the result of homoplasy. Furthermore, because antipredator features may be secondarily lost, nettle, and gelatine caterpillars may be the result of multiple gains or losses of spines.
Hence, our objectives were as follows: (a) to reconstruct a well-supported phylogeny of the Limacodidae using a multigene dataset and (b) to trace the evolution of spines by optimizing character states of slug caterpillars with and without spines on this phylogenetic framework. We also comment on the systematic relationships of the Limacodidae. Most of the taxa included in this study were reared from samples collected from Asia, but we also include material from Australia, North America, and South America.

| Molecular data
Total genomic DNA was extracted from 1 to 3 legs of each specimen using a commercial DNA extraction kit (Gentra Puregene Tissue kit, Qiagen) following the manufacturer's protocol. The polymerase chain reaction (PCR) was used to amplify the following five gene fragments: cytochrome oxidase subunit I (COI), D2 region of the 28S ribosomal sequence, 18S ribosomal sequence, elongation factor-1 alpha (EF-1α), and partial sequences of the wingless gene. The first mentioned fragment is encoded in the mitochondrial genome, whereas the remaining four markers are part of the nuclear genome. These genetic markers are phylogenetically informative and commonly used for resolving the systematics of the Lepidoptera (Chalwatzis, Baur, Stetzer, Kinzelbach, & Zimmermann, 1995;Lee & Brown, 2008;Lo et al., 2015;Mutanen, Wahlberg, & Kaila, 2010;Niehuis et al., 2006;Regier et al., 2013Regier et al., , 2009Simon et al., 1994;Wahlberg & Wheat, 2008;Zaspel et al., 2016). A list of primers used for generating sequence data from the targeted loci is given in
TA B L E 1 List of species used in the phylogenetic analysis for this study, their broad geographical distribution, larval character states A-D (A = spines present after second instar; B = spines present after second instar but reduced in late instars; C = spines and setae absent in all instars; D = spines absent but numerous setae present after second instar), and GenBank accession numbers Altha melanopsis

| Larval morphology
We collected eggs and larvae for most species to record larval character states. Some eggs were obtained from females collected from light traps, while other eggs and larvae were collected directly from the field. Eggs and larvae were brought back to the laboratory and assigned rearing records, adopting the system used by Powell and De Benedictis (1995
Based on previous studies (Battisti et al., 2011;Epstein, 1996 State B: Spines present after the second instar (Figure 1j,n,r,u), but almost all spines are lost or reduced in late instars (Figure 1k,o,s,v); when the spines are reduced, they are tiny and vestigial ( Figure   1v). A few setae are present on pairs of protuberances on each segment in the first instar (Figure 1i,m,q).
State C: Spines absent in all instars (Figure 1t,w,x). Further, the setae in the first instar are also vestigial, such as Belippa horrida (Epstein, 1996).
State D: Spines absent; numerous setae are present on tubercles, which can be pulled off after the second instar ( Figure 1p); a few setae are present on pairs of protuberances on each segment in the first instar (Figure 1l).

| Character evolution analyses
The character evolution of larval spine variation was reconstructed on the maximum clade credibility tree using the Mk1 evolutionary model as implemented in Mesquite (version 3.2) (Maddison & Maddison, 2017).

| Phylogenetic patterns
The

| Character evolution of spines
The evolutionary reconstruction of spines in limacodid caterpillars indicated that the ancestral state was most likely larvae with spines present from second instar to final instar (character state A) (Figure 3, Node 1: TA B L E 2 List of primers used for generating sequence data for the five genetic markers

| D ISCUSS I ON
Our molecular study provides a robust phylogeny of the Limacodidae.
The well-supported phylogenetic framework allows us to reliably reconstruct the character evolution of spines throughout the entire F I G U R E 2 Phylogenetic trees of the Limacodidae based on the combined dataset constructed with: (a) partitioned Bayesian Inference; (b) partitioned Maximum Likelihood using the GTR + Γ+I substitution model. Branch lengths are proportional to inferred nucleotide substitutions, with values above nodes representing posterior probabilities (a) and ML bootstraps (b). Optimal topologies recovered by BI and ML were congruent. Six major lineages were recovered, which are indicated by different colors. Zoogeographic regions are represented in different colors on terminals, as per legend larval stage, to test previous hypotheses regarding the evolution of slug caterpillars, and to infer the potential mechanisms of homoplasy in limacodids.

| Character evolution and morphological homoplasy
According to the phylogeny reconstructed in this study, limacodids evolved from a common ancestor in which the larval type is consistent with Cock's (1987) hypothesis that nonstinging types of slug caterpillars evolved from nettle caterpillars. Although the pattern contrasts with the larval character evolution of Zaspel et al. (2016), it must be emphasized that branch support for many of the basal nodes in that phylogenetic study was low and hence ancestral reconstructions were at best preliminary.
Spines in the Limacodidae are considered to be an adaptive response to predation (Murphy et al., 2010). Our phylogeny indicates that this defense strategy evolved early in the origin of the family, and the trait is widespread across lineages 1 and 4-6 ( Figure 3).
Therefore, the independent losses of poisonous spines (homoplasy) raise the interesting question as to why have some larvae evolutionary lost their toxic antipredator mechanism? Gelatine caterpillars avoid predation through crypsis or masquerade, but it remains to be determined what mechanism may have driven this type of defense strategy. Here, we propose several potential mechanisms (hypotheses) for spine reduction in slug caterpillars.
The first hypothesis is that spines get lost or reduced because they confer no advantage below a certain size threshold. It has been demonstrated that defensive characters such as warning coloration are more effective when displayed in insects with large bodies (Forsman & Merilaita, 1999;Hossie, Skelhorn, Breinholt, Kawahara, & Sherratt, 2015). For example, defensive eyespots are effective in big caterpillars, but costly in small caterpillars, because they enhance detectability without providing a protective advantage in small caterpillars (Hossie et al., 2015). In tree-feeding insects, avian predation risk increased with larger prey body size (Remmel, Davison, & Tammaru, 2011;Remmel & Tammaru, 2009). Therefore, slug caterpillars with small body size (e.g., Quasinarosa corusca) may be hard to detect, so that the cost of producing spines and toxins may be higher than the benefit of avoiding predation in smaller taxa.
The second hypothesis is that there has been a change in predator pressure. Predators (e.g., insectivorous birds) eat aposematic prey in a selective manner according to their levels of hunger and the presence of alternative prey (Cott, 1940;Ruxton et al., 2004). When limacodids expand their range or enter new adaptive zones, such as in low diversity biomes (e.g., high mountain or desert habitats), with potentially higher levels of predator pressure and less alternative prey, nettle caterpillars may be too obvious to survive and cryptic larvae without spines may be selected for.
The third hypothesis is that slug caterpillars without spines may be physiologically more suited to dry environments, such as deserts, seasonal savannas, and alpine woodlands (Leuschner, 2000).
According to previous studies (Battisti et al., 2011;Cock et al., 1987;Epstein et al., 1999;Hossler, 2010;Kano, 1977;Kawamoto & Kumada, 1984), spines on nettle caterpillars consist of multiple cells, and spines are usually arranged on tubercles. Slug caterpillars with spines on tubercles have higher surface area to volume ratios than slug caterpillars without spines and tubercles. Surface area to volume ratios may influence water balance in ectotherms (Ashton, 2002;Bidau & Marti, 2008). For example, the tropical rain frog, Eleutherodactylus coqui, reduces water loss by adjusting posture and activity to control the exposed surface area (Pough, Taigen, Stewart, & Brussard, 1983;Vitt & Caldwell, 2013). By analogy, slug caterpillars without spines with lower surface area to volume ratios may be more suited to dry environments. In a previous study, it has been observed that nettle caterpillars are distributed more in tropical areas and gelatine caterpillars are distributed more in temperate areas (Zaspel et al., 2016).
In addition to adaptation to similar local environments, because genetic or developmental constraints limit the generation of phenotypic variations (Brakefield, 2006;Hall, 2007;Wake et al., 2011), the reappearance of similar features in organisms may result from different selective pressures (Hall, 2007). For example, pelvic reduction in stickleback populations, which are sympatric with various fish and bird predators, may be triggered by low calcium ion concentration (Giles, 1983); in Paxton Lake with a high calcium ion level and in some Alaskan Lakes with lack of native predatory fishes, stickleback populations have similar pelvic vestiges (Bell et al., 1985;Larson, 1976). Therefore, homoplasy of pelvic reduction in sticklebacks is more likely to be caused by different selective pressures, low calcium ion concentration and lack of native predatory fishes, in different lakes (Bell, 1987). Furthermore, homoplasy is common with reduced characters especially for complex characters, which may have low probability of origin but can be lost or reduced by the action of a few genes (Culver & Pipan, 2016;Cunningham, Omland, & Oakley, 1998;Maddison, 1994;Sackton et al., 2019). In this study,

| Systematic considerations
In the inferred phylogenetic tree of the Limacodidae, we identified six lineages (Figure 2 character state A. The structure of the spine in lineage 1 is the same as that in lineages 4-6, which is formed by trichogen cells that line up with the epidermal cells (Kawamoto & Kumada, 1984), although the numbers of spines on each segment (Figure 1d,h) are fewer than those in lineages 4-6 ( Figure 1b,c,f,g).
In lineage 2, three taxa comprise a monophyletic group that is characterized by hairy monkey slug caterpillars (character state D). The clade includes Isochaetes sp. and Phrixolepia inouei, which emerged as sister taxa. The geographical distribution of Isochaetes is in eastern North America, Central America, and northern South America, whereas the distribution of Phrixolepia is mainly in eastern Asia (Ratnasingham & Hebert, 2007). The disjunction between North America and eastern Asia has been reported for many animal and plant taxa (Espeland et al.2015;Nordlander, Liu, & Ronquist, 1996;Peña, Nylin, Freitas, & Wahlberg, 2010;Tiffney, 1985;Wen, 1999 With the exception of Chalcocelis albiguttatus, all other taxa from Australia (seven species representing six genera) comprised a monophyletic group within lineage 6 ( Figure 3: Node 5). Although the clade was not strongly supported, it may be improved by greater taxon sampling of the fauna of the continent. The topology and relative branch lengths indicate that most limacodids in Australia evolved relatively recently. Moreover, the Australian lineage is nested within a set of predominantly Asia lineages (lineages 4-6), which suggests that the origin of these limacodids is not in Australia. Further taxon sampling of the family and divergence times using a molecular clock are needed to estimate deeper biogeographic patterns to test this hypothesis.

ACK N OWLED G M ENTS
We appreciate David Wagner's and an anonymous reviewer for their careful review and insightful suggestions. We thank Shou-

CO N FLI C T O F I NTE R E S T
None declared.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in GenBank at https ://www.ncbi.nlm.nih.gov/genba nk/, accession numbers in Table 1.