­Comparative spigot ontogeny across the spider tree of life

Spider husbandry

Rearing conditions and lab colony maintenance for Hogna carolinensis and Dolomedes tenebrosus follow Mallis & Miller (2017). The founding female H. carolinensis was collected in Bernalillo County, New Mexico (D. Lightfoot, 10-Sept-2014) carrying second instar spiderlings on her abdomen. A gravid D. tenebrosus with an egg sac was sent to us from Bedford County, Virginia (K. Benson, 29-June-2014). When possible, two to three replicates of each instar were randomly sacrificed for scanning electron microscopy (SEM) imaging to account for variation among individuals of the same instar. After developing through the first few instars of both species, colony survival strongly declined. Subsequently, single samples were collected at each instar for D. tenebrosus, while after the seventh instar in H. carolinensis, a single female was followed to adulthood (twelfth instar). The founding females, as well as instar vouchers, were deposited at the Museum of Southwestern Biology, Division of Arthropods (MSBA 50049–50070).

SEM preparation, imaging and spigot mapping

We dissected in 100% EtOH, and then critical point dried using a Denton DCP-1 critical point dryer, and mounted specimens of each instar for D. tenebrosus and H. carolinensis onto standard SEM stubs using a combination of copper tape and non-conductive glue at California Academy of Sciences (CAS). Finally, all SEM stubs were coated in gold/palladium using a Cressington Sputter coater 108 (6002, 6006 series) that used Argon gas to facilitate coating. At CAS, we obtained SEM images on the Hitachi SU-3500 scanning electron microscope. Up to 20 views of each instar spinning field for both species were captured, covering all spinnerets for two replicates of each instar. These SEMs were used to create spigot maps which were translated into a spigot ontogeny dataset for each species (methods outlined in Coddington (1989) and Griswold et al. (2005); see Table 1). Spigot mapping allows for notation of type, functionality and placement of silk spigots on each spinneret, as well as tracking the growth of the spinning fields from instar to instar and the appearance of interesting spigot formations.

Spigot ontogeny datasets

We compiled a large spigot ontogeny dataset of 22 species comprising thirteen spider families using previously published studies and unpublished datasets (Table 2; Wąsowska, 1977; Yu & Coddington, 1990; Hajer, 1991; Townley & Tillinghast, 2009; Dolejš et al., 2014; R Carlson & CE Griswold, 1996, unpublished data). Three of those sources came from our lab colonies not only for D. tenebrosus and H. carolinensis, as reported here but also the cribellate zoropsid, Tengella perfuga (RE Alfaro, 2017, unpublished data). The dataset included the appearance, type, and number of specific spigots on each spinneret (Fig. 1, Table S1). While some studies reported tartipores (remnant structures of spigots from previous instar spigots) and others the presence of nubbins (non-functioning spigots), all studies reported the number of the seven common (shared) spigot types across all the species of this study: 1. MAP gland spigots, 2. PI gland spigots, 3. mAP gland spigots, 4. AC gland spigots on the PMS, 5. CY gland spigots on the PMS, 6. AC gland spigots on the PLS and 7. CY gland spigots on the PLS. We standardized the final data set to include these ‘standard 7’ spigots (Table S1).

Variables of interest

From the standardized datasets, we selected five independent or predictor variables for phylogenetic comparative analyses, and chose twelve dependent or response variables for the final analyses. The independent variables were (Table 2):

1. Strategy: Foraging strategy scored as 0: web builder or 1: non-web builder.

2. Specific: Specific foraging strategies, scored as 1: sit & wait, 1.5: ambush, 2: sit & pursue, 2.5: stalking, 3: active, 4: sheet web, 4.5: funnel web or 5: orb web. Sit & pursue, while similar to sit & wait, is defined as an intermediate state between sit & wait and active hunting (Schmitz & Suttle, 2001; Miller, Ament & Schmitz, 2014). With a sit & wait strategy, the otherwise motionless spider grabs the prey when it comes within striking distance, whereas with sit & pursue, the otherwise motionless spider runs after the prey when it comes within several centimeters proximity and pursues the prey until captured (Schmitz & Suttle, 2001; Miller, Ament & Schmitz, 2014).

3. Silk: Main type of silk used, scored as 1: none/MAP dragline, 1.5: burrow, 2: aciniform sheet, 3: cribellate, or 4: viscous silk, i.e., that produced from the aggregate gland spigots.

4. Type: A measure of the variety of spigots the species possessed beyond what we called the standard 7. These scored as 1: standard 7, 1.5: standard 7 + modified piriform gland spigots on the ALS, 2: standard 7 + modified spigot on PLS, 2.5: standard 7 + modified spigot on PLS with two ‘flankers’ + cribellum, 3: standard 7 + aggregate and flagelliform gland spigots on the PLS or 4: standard 7 + cribellum, paracribellar gland spigots on the PMS and a pseudoflagelliform gland spigot on the PLS.

5. Instar: Maximum number of instars that a species goes through to reach adulthood (females).

The twelve dependent variables were continuous and were the average number of spigots for each of the ‘standard 7’ found in all the spider species spinning fields. We focused on the adult female instar as well as the second instar, i.e., when all spiders emerge from the egg sac. Specifically, for adult females, the seven dependent variables were average number of spigots for: 1. MAP gland spigots on the ALS, 2. PI gland spigots on the ALS, 3. mAP gland spigots on the PMS, 4. AC gland spigots on the PMS, 5. CY gland spigots on the PMS, 6. AC gland spigots on the PLS and 7. CY gland spigots on the PLS. For second instars, the five dependent variables were the average number of spigots for: 8. MAP gland spigots on the ALS, 9. PI gland spigots on the ALS, 10. mAP gland spigots on the PMS, 11. AC gland spigots on the PMS and 12. AC gland spigots on the PLS. Cylindrical gland spigots are only found in adult female spiders, which is why they are not included with the second instars.

Phylogenetic comparative analyses have typically dealt with at least two continuous variables. Methods have improved to accommodate the increase in Type 1 Error associated with discrete variables such that analyses like phylogenetic generalized least square models (PGLS) can be robustly performed (Graber, 2013; Maddison & FitzJohn, 2015). While theoretically repeated-measures or a factorial ANOVA are possible in a phylogenetic context, the methods of incorporating them into phylogenetic comparative analyses are not yet developed (Guo et al., 2007). Thus, our time series study could not be analyzed as a whole unit. To gain a picture of the potential effects of ontogeny, we performed several PGLS analyses on the adult female spigot numbers as well as those of the second instars. These include the stage at which all spiders emerge from the egg sac and are the most similar in condition (second instars) to the stage when the most diversification and growth in spigots appears (adult female stage).

Phylogeny

We used the topology of the recently published Araneae Tree of Life (AToL) (Wheeler et al., 2016) for the 22 species included in our study. At the time of publication of this manuscript, the AToL sequences have not been publicly released. Thus, while we were not able to generate a time-calibrated tree or one with branch lengths equal to a rate of molecular evolution, we could use the published topology and create two phylograms. One phylogram had all branch lengths equal to 1 and the other was an ultrametric tree with the same topology (Fig. 2), using the ape (version 4.1) package in R (Paradis, Claude & Strimmer, 2004). Both trees were used in model testing for phylogenetic comparative analyses and ultimately, we used the ultrametric topology in our final analyses.

Phylogenetic comparative analyses

All analyses were performed in R using RStudio (version 3.2, R Core Team, 2016; RStudio, 2017). First, we tested both sets of dependent and independent variables for phylogenetic signal using the phytools (version 0.5–64) package in R (Revell, 2012). All independent variables had strong phylogenetic signal with Blomberg’s K (Blomberg, GarlandJr & Ives, 2003) and Pagel’s λ (Pagel, 1999) being close to 1 and with significant p-values while most of the dependent variables did not show phylogenetic signal (Table 3). Despite this, given the strong signal in the predictor variables, we decided to proceed with PGLS analyses. Because the K values were so close to 1, a Brownian motion model of evolution was the best fit. We further tested various models of evolutionary rate, and confirmed the Brownian model of evolution being the best fit (Results S1). Next, using both trees, we tested PGLS models using generalized least squares method of model selection with single term up to all five independent variables included. We used the R packages, nlme and MuMIn (version 1.15.6) to test model fit (Bartón, 2016; Pinheiro et al., 2017). PGLS analyses required the use of the ape, geiger (version 2.0.6), nlme (version 3.1–131) and phytools packages in R (Paradis, Claude & Strimmer, 2004; Harmon et al., 2008; Revell, 2012; Pinheiro et al., 2017). With delta AICc values equal to zero or very close to zero between terms, and a significant p-value associated with the single term Instar, we determined that single term models were the best fit with the trees and the datasets for PGLS as there was no significant effect of adding additional terms (Table 4, Results S2). We ran a total of 35 PGLS analyses covering each of the adult female dependent variables (average numbers of specific silk spigots) and 25 PGLS analyses covering each of the second instar dependent variables using the ape and geiger packages in R (Paradis, Claude & Strimmer, 2004; Harmon et al., 2008). We did this with both trees, and while the specific analyses results differed, the main conclusions did not. We also performed ANOVA analyses on the independent variable means derived through the PGLS models. Because all the independent variables were discrete the ANOVA analyses served to corroborate the PGLS results. Here we report the significant PGLS coefficient of correlation results (Table 5) using a maximum likelihood approach, with a Brownian model of evolution and the ultrametric tree. The 60 full PGLS analyses and ANOVA results using the ultrametric tree are available in Results S3 and S4.

Ancestral character estimation

Finally, we used the ultrametric tree and the ace function in the ape package in R (Paradis, Claude & Strimmer, 2004) to conduct ancestral character estimation on the diverse, singular spigots found on the PLS to explore the unresolved question of whether these spigots are homologous structures or not. The spigots of interest were the modified spigot (MS) found in some cribellate and ecribellate spiders, the pseudoflagelliform gland spigot (PF) found in cribellate orb weavers and the flagelliform gland spigot (FL) found in viscous orb weavers. We added five additional taxa, using data derived from adult female SEM images in Griswold et al. (2005) to allow for broader taxon sampling deeper into the phylogeny. These species were the cribellate Hypochilus pococki N.I. Platnick (1987; Hypochilidae), Kukulcania hibernalis (N.M. Hentz, 1842; Filistatidae), Thaida peculiaris F. Karsch (1880; Austrochilidae), Megadictyna thilenii F. Dahl (1906; Megadictynidae) and ecribellate Nicodamus mainae M.S. Harvey (1995; Nicodamidae). The Nicodamoidea (Megadictynidae plus Nicodamidae) are sister to all Araneoidea, while the other families are sister to the Araneoidea + RTA clades (Wheeler et al., 2016). We used a maximum likelihood method with a model of the weighted rate matrix of substitutions for these spigots (Table 6). As we describe below, we constrained our rate matrix (Table 6) that we used to model substitution rates across branch lengths based on prior knowledge about historical possession of the flagelliform gland, pseudoflagelliform gland and modified spigots in each taxon (Table 6) (Platnick & Griswold, 1991; Griswold et al., 2005; Dimitrov et al., 2016; Wheeler et al., 2016). Some of the more recently derived clades within the Araneoids (Dimitrov et al., 2016) have lost flagelliform gland spigots (i.e., Mimetidae, Arkyidae). We ranked the transition to Loss (of spigots) as 1 substitution and that to none, pseudoflagelliform, or modified spigots (not observed anywhere in this lineage) as 0 substitutions, or no likelihood (Table 6). Because they are found in cribellate groups ancestral to both the RTA and Araneoid clades, we ranked all transitions for modified spigots to the other four states as 1 substitution (Table 6). One interpretation of Dollo’s Law is that it is easier to lose a structure than to re-evolve it (Dollo, 1893) and this influenced how we weighted the remaining matrix for Pseudoflagelliform, Loss, and the possible None (no spigot) cases. Since one of the cribellate ancestral species, Kukulcania hibernalis, does not possess anything corresponding to a modified or pseudoflagelliform gland spigot, and is classified as ‘None’, we conservatively allowed for a single substitution from None to Flagelliform, Modified and Pseudoflagelliform in our rate matrix (Table 6). Finally, since the cribellate orb weaving sister group to the RTA clade possesses a pseudoflagelliform gland spigot (PF) and RTA clade members possess modified (MS) or no spigots, we allowed for 0 substitutions in this direction. We also allowed for a substitution rate of 1 between Pseudoflagelliform and Flagelliform as we do not know whether the ancestral orb weaver possessed pseudoflagelliform or modified spigots (Bond et al., 2014; Garrison et al., 2016; Wheeler et al., 2016). We then plotted the likelihoods of states at each node on the ultrametric phylogeny (Fig. 3, Results S5).

Additionally, we performed ancestral character estimation analyses for four of the five independent variables from the PGLS analyses on the ultrametric tree with increased taxa: Strategy, Specific, Silk and Type (Results S5). The constrained rate matrices for each variable along with explanations for the given matrices are provided in Table S2. We did not include maximum instar, as those data were not available in the literature for the five additional taxa. Also, due to its variability within a species (see Table S1), ancestral character estimation may not be a useful tool to use to understand maximum instar in silk evolution.

Spigot ontogeny of Dolomedes tenebrosus and Hogna carolinensis

All instars of D. tenebrosus (instars 2–13) and more than half of the instars (instars 2–7 and 12) of H. carolinensis were observed and sampled from lab colonies for SEM imaging to assess the spigot ontogeny of the full spinning field of these lycosoid spiders (Homann, 1971; Griswold, 1993; Polotow, Carmichael & Griswold, 2015). Dolomedes tenebrosus reached adulthood in thirteen instars while H. carolinensis reached adulthood in twelve (Table 1).

Anterior lateral spinnerets

Both species possessed two MAP gland spigots, except for the adult male stage (Table 1, Fig. 4). Piriform gland spigots increased in number to adulthood (Fig. 4). However, in adult male D. tenebrosus, PI gland spigots decreased from the penultimate instar, which led to sexual dimorphism (Table 1). It was not until later instars of D. tenebrosus and H. carolinensis that the number of PI gland spigots increased in greater magnitude from instar to instar (Figs. 4A–4E, Table 1). We also observed sensilla (sensory pores) in the MAP fields of both species (Fig. 5).

Posterior median spinnerets

Both species possessed two mAP gland spigots, except for the adult male stage of D. tenebrosus (Table 1, Fig. 6). CY gland spigots did not appear until the adult female and probably penultimate female instars. Penultimate D. tenebrosus possessed at least one pre-cylindrical gland spigot (Fig. 6F). D. tenebrosus adult females bore many more CY gland spigots than H. carolinensis (Fig. 6C). In D. tenebrosus, AC gland spigot numbers did not change at 4 AC gland spigots for instars 2/3, then at 5 AC gland spigots through instars 4/5/6, and then 8 AC gland spigots for instars 8/9/10 (Table 1, Figs. 6D, 6E). Both males and females lost AC gland spigots in the final molt to adulthood. In H. carolinensis, AC gland spigots dropped in number at instar 4 and remained at 3 AC gland spigots for the next two instars (Table 1, Figs. 6A, 6B). By the adult instar, the number was far greater in the female H. carolinensis than the female D. tenebrosus (Table 1, Fig. 6C).

Posterior lateral spinnerets

Cylindrical gland spigots also appeared in the penultimate female stage with the three pre-spigots visible (Fig. 7F). Adult female H. carolinensis possessed a single CY gland spigot compared to the 28 CY gland spigots in adult female D. tenebrosus (Table 1). Aciniform gland spigots presented different trends between the two species (Table 1, Figs. 7 and 8). In D. tenebrosus, AC gland spigot numbers slowly increased with the same number persisting for 2–3 instars, then increasing (Fig. 7). We also observed loss of spigots in the final molts to adulthood in both males and females (Table 1). In H. carolinensis, a sharp decrease to 3 AC gland spigots occurred in instar 4 and persisted through instar 5, then increased to 7 AC gland spigots in instar 6 persisting through instar 7 (Table 1, Fig. 8). The female H. carolinensis possessed more AC gland spigots than the female D. tenebrosus. In both species, a larger spigot was tentatively identified as a ‘modified spigot’ (see Griswold et al., 2005: 61; character 96), with a potential pre-modified spigot observed in the penultimate female stage of D. tenebrosus. These made no other appearance in the ontogeny the spinning fields of either species (Table 1, Figs. 7 and 8).

Phylogenetic comparative analyses

We compiled a full spigot ontogeny dataset of 22 species representing thirteen spider families and scored the five independent, predictor variables: Strategy, Specific, Silk, Type, Instar (Table S1, Table 2). After determining that the independent variables had strong phylogenetic signal, whereas only three of the twelve dependent response variables did, we decided to proceed with the phylogenetic generalized least squares analyses (Table 3). After thorough tree and model selection analyses, we determined that single term models were the best fit along with the ultrametric tree and a Brownian model of evolution (Table 4, Fig. 2). We ran a total of 60 PGLS analyses (Results S3, S4) and here report the significant results of those analyses (Table 5).

In the adult female analyses, we found that Instar was a significant predictor in a few cases. For piriform gland spigots (ALS), the coefficient of correlation between Instar and average number of spigots was significant (Table 5). This was also the case for aciniform gland spigots of the PMS (Table 5). For female PI gland spigots, Type, or the variety of spigots possessed, was also a significant predictor of Average number of PI (Table 5). Of interest, the female mAP gland (PMS) spigots did show strong phylogenetic signal and Strategy was a significant predictor of Average number of mAP gland spigots. This means that, after the bias of the correlation due to phylogeny was accounted for through the PGLS, the species classification of either being a web builder or not had a strong coefficient of correlation with the average number of mAP gland spigots possessed (Table 5).

In the second instar analyses, we found that Instar was also a significant predictor or independent variable. For MAP gland (ALS) spigots in second instars, the coefficient of correlation between Instar and Average number of MAP gland spigots was significant (Table 5). Spiders possessed one, two or no MAP gland spigots at the emergence from the egg sac (second instar) (Table S1). For the mAP gland (PMS) spigots, which did show strong phylogenetic signal in preliminary testing, Instar was a significant predictor of the Average number of mAP gland spigots (Table 5). Second instar spiders possessed two, one or no mAP gland spigots at this first stage outside of the egg sac (Table S1).

Ancestral character estimation of singular spigots on the PLS

We performed a maximum likelihood ancestral character estimation (ACE) of potentially homologous spigots producing axial lines, i.e., the singular spigot on the PLS. We used the ultrametric tree and overall scoring of whether a species possesses a spigot and if so, which spigot it is: Flagelliform, Modified, None (never historically possessing a singular spigot), Pseudoflagelliform, or Loss (previously possessing a singular spigot) of a spigot. We used a constrained rate matrix for the five states as our model of substitution rates for the ACE analyses (Table 6). The log likelihood value for the analysis was −46.023. Within this analysis we calculated the scaled likelihoods of states at the root (Flagelliform gland spigot: 0.000, Modified spigot: 0.365, None (no singular spigot that produces axial lines): 0.635, Pseudoflagelliform gland spigot: 0.000, Loss: 0.000), as well as the other nodes of the phylogram produced (Fig. 3, Results S5). In our color-coded results, the predominately red clade was the Araneoidea, which includes viscous orb weavers (FL spigot) and the pirate spiders (Mimetidae) which do not possess a FL spigot as adults (Loss) (Fig. 3). However, as very young juveniles mimetids possessed vestigial FL and aggregate gland spigots on the PLS (Table S1, Townley & Tillinghast, 2009). The predominately green clade (MS spigot) was the RTA clade, with loss of MS spigot (black) in many of the tip species, while just sister to that were the yellow clades (Pseudoflagelliform), which included the cribellate sheet (Phyxelida) and cribellate orb weavers (Hyptiotes). Additional plesiomorphic cribellate taxa, possessing MS spigots (green, Thaida, Hypochilus) or no unique spigot on the PLS (blue, Kukulkania) were included for deeper phylogenetic sampling (Fig. 3). Finally, the results suggested that, at the ancestral root node, the PLS singular spigot was more likely to have been a modified spigot than a pseudoflagelliform gland spigot, and that the ancestor either possessed a modified spigot or no singular spigot (Fig. 3, Table 4).

We followed a similar procedure for maximum likelihood ancestral character state estimation for four of the five independent variables: Strategy, Specific, Silk and Type. The constrained rate matrices used as the model for the analyses are given in Table S2. The resulting state likelihood values for each node are also given in Results S5 with the previous ACE node results. For the variable Strategy (foraging strategy), the ancestral character estimation model has a log likelihood of −14.73. At the root, the ancestral state was “Web” or web building (likelihood = 1.0) rather than no web or loss of web (Fig. 9, Results S5). Therefore, the ancestral state for the araneoids and lycosoids was a web builder. Specific (specific foraging strategy) ACE model has a log likelihood of −58.65. At the ancestral node, the states are for web building (Fig. 10). The ancestor of lycosoids and araneoids likely spun a funnel web (likelihood = 0.812) or a sheet web (likelihood = 0.084); an orb web (likelihood = 0.105) is less likely (Fig. 10, Results S5). The ACE model for Silk (main type of silk expressed) has a log likelihood of −44.61. The ancestor of araneoids and lycosoids produced cribellate silk (likelihood = 1.0; Fig. 11, Results S5). Finally, the ACE model for Type (variety of silk spigot types) had a log likelihood of −36.14. The ancestral state at the root node was undetermined, but the two states with the highest likelihoods were cribellate spiders possessing either Standard 7 + Cribellum + Paracribellar spigots (likelihood = 0.398) or Standard 7 + Cribellum + Pseudoflagelliform + Paracribellar spigots (likelihood = 0.433). All additional five states had much lower likelihoods (Fig. 12, Results S5).

Spigot ontogeny of D. tenebrosus and H. carolinensis

Both D. tenebrosus and H. carolinensis are large-bodied lycosoids. Dolomedes tenebrosus belongs to the Pisauridae family, the fishing or nursery web spiders, and employs a sit & wait foraging strategy (Table 2) whereas H. carolinensis, belonging to the diverse Lycosidae family, employs a sit & pursue strategy (Table 2) (Schmitz & Suttle, 2001; Miller, Ament & Schmitz, 2014). While both have a similar number of instars to adulthood, the two species differ dramatically from each other in the loss and regain of aciniform gland spigots on the PMS and PLS (Table 1, Figs. 6, 7 and 8). In D. tenebrosus, AC gland spigot numbers on the PMS do not change for multiple instars slowly increasing after two to three instars with a dramatic increase in the final molts. In H. carolinensis, AC gland spigot numbers drop and remain low for two to three instars before gradually increasing, until a similar dramatic increase in number with the final molt (Table 1, Fig. 6). The same interesting trend was also observed for aciniform gland spigots on the PLS in both species (Table 1, Figs. 7 and 8). These trends were not observed in other lycosoids previously studied including other members of the Lycosidae (Wąsowska, 1977; Dolejš et al., 2014). H. carolinensis do not form webs. D. tenebrosus instars were observed in the lab to build silk scaffolding where they rested and at times fed (RE Alfaro, 2017, unpublished data). It is possible that this difference in silk use and foraging strategy between H. carolinensis and D. tenebrosus, especially the lack of regular web building, could account for the trends we observed with AC gland spigot numbers. Both species were observed on a few occasions to wrap prey items after a preliminary bite before returning to bite again. This is an ancient behavior seen in many other taxa including the Araneoidea, Phyxelididae, and other members of the RTA clade. It is also possible that they have evolved alternative uses for aciniform silk or do not need to produce this silk until the adult instars when numbers of spigots on both spinnerets increase.

Sexual dimorphism was observed in PI gland and AC gland spigots in D. tenebrosus and this has also been recorded in other lycosoids and araneid spiders (Wąsowska, 1977; Townley & Tillinghast, 2009; Dolejš et al., 2014). We were not able to rear our H. carolinensis males to adulthood in the lab colony, but given the trends in other Lycosidae, we would expect to observe sexual dimorphism as well (Dolejš et al., 2014). In most spiders, regardless of lineage or cribellate or ecribellate status, adult males lose (abort) spigots of all types in the final molt: this is likely due to the shift in life history strategy of abandoning webs or territories to actively forage and look for females (Wąsowska, 1977; Yu & Coddington, 1990; Hajer, 1991; Townley & Tillinghast, 2009; Dolejš et al., 2014; R Carlson & CE Griswold, 1996, unpublished data; RE Alfaro, 2017, unpublished data).

Adult females of both species possessed a potential modified spigot (MS) on the PLS. A modified spigot was not reported for other lycosoids in previous studies (Wąsowska, 1977; Dolejš et al., 2014). Modified spigots have been reported as singular or, in some cases of cribellate silk users among members of the RTA clade, which includes lycosoids, the MS may occur with flanking spigots (Griswold et al., 2005; RE Alfaro, 2017, unpublished data). It is possible that this modified spigot in lycosoids could be homologous to the pseudoflagelliform gland spigot (PF) observed in cribellate lineages sister to the RTA clade, e.g., Phyxelida tanganensis (E. Simon & L. Fage, 1922) and Hyptiotes paradoxus (C.L. Koch, 1834) (Table S1; Peters, 1984; Peters, 1995; Griswold et al., 2005; Eberhard, 2010; Eberhard & Barrantes, 2015; R Carlson & CE Griswold, 1996, unpublished data). In contrast, it should be noted that in species with modified spigots, adult males typically possess an MS nubbin, which we did not observe in male D. tenebrosus (Ramírez, 2014, Table 5). To confirm these observations in D. tenebrosus and H. carolinensis, replicate adult female and male specimens are necessary.

Phylogenetic comparative analyses

Six analyses had significant correlation coefficients after correcting for the bias of shared evolutionary history, suggesting that Instar, Strategy, and Type (spigot variety) are good predictors of spigot number in spiders (Table 5). Although most of the PGLS analyses found no significant correlation of the five independent variables with the 12 dependent variables of spigot numbers, this is not altogether surprising considering the PGLS analyses remove the bias of correlation due to shared evolutionary history. It is possible that with broader taxon sampling deeper in the phylogeny, trends may emerge beyond those explained by phylogenetic signal. As our analyses currently cannot include as a variable the full ontogeny dataset per species, it is also more likely with future developments of more complex statistical analyses within a phylogenetic context that will allow us to include the full ontogeny picture of each species, we will be able to gain a better understanding of what is driving silk spigot evolution in spiders.

Maximum number of instars (Instar) served as a proxy for body size or body condition in each spider species. Within a species with variation in maximum instar, we observed higher numbers of spigots in the individuals of older instars (Table S1). Finding a significant correlation between Instar and Average number of spigots in adult female aciniform and piriform gland spigots is not unexpected, considering these spigots increase in number with each instar (Table S1, Table 5). The more nutrition a juvenile spider consumes in one instar influences how much growth occurs in the molt to the next instar. Spiders with a steady food supply may invest in an increased number of instars to ensure better body condition at the adult stage, which could lead to a trend of increase in spigot numbers (RE Alfaro, 2017, unpublished data). Piriform gland spigot numbers were also positively correlated with Type (variety of spigots possessed) (Table 5). Piriform gland silk is used as a cement for other silk fibers, particularly the major ampullate gland fibers in dragline silk (Garb, 2013). In orb weavers, piriform gland silk is used to cement the structural lines together and to the substrate whereas in wandering spiders it cements the dragline to the substrate to prevent the spider from falling (Garb, 2013). Orb weavers and lycosoids possess different types of silk spigots and have different uses for the shared spigots they possess (Table S1, Table 2); for example, MAP gland silk may be used to construct aerial frames (orb weavers) vs. surface dragline (lycosoids); aciniform gland silk may be used as web material in wolf spiders who spin funnel webs (the genus Hipassa) rather than as a prey wrapping material (orb weavers) (Mathew, Sudhikumar & Sebastian, 2011; Garb, 2013). It is possible that these differences in use are due to the types of spigots they possess and their differing foraging strategies (web building vs. predominately active hunting) and are what is causing this positive correlation we observed (Table 5). Finally, in adult female PGLS analyses, Strategy (web vs. non-web) was a significant predictor of the number of mAP gland spigots on the PMS. The coefficient of correlation was small, but when we look at the full ontogeny data (Table S1), we see a clear differentiation between araneoids and the others. In the adult female stage, araneoids lose one mAP gland spigot and retain one functional spigot, whereas in the other groups, from the lycosoids to the cribellate web builders, all female spiders retain the two mAP gland spigots that they possessed throughout their ontogeny (Table S1). The clear correlation between strategy as a predictor and mAP number as a response is expected since all araneoids possess one mAP spigot and the remaining spider groups possess two (Correlation coefficient: 0.291, t = 2.448, p = 0.024: Table 5). This may have withstood the phylogenetic correction of the PGLS, because several of the non-araneoids within our study also spin webs.

The PGLS analyses for second instars were largely non-significant. The lack of significant correlation between predictor and response variables may be due to second instars being more similar. Many spiders start out with the same general number of spigots upon emergence from the egg sac (second instar) and differentiation between foraging strategies may not be apparent at this instar, e.g., some second instar web builders do not spin webs (Table S1, Hajer, 1991; Barrantes & Madrigal-Brenes, 2008; RE Alfaro, 2017, unpublished data). However, both MAP gland spigots on the ALS and mAP gland spigots on the PMS were significantly correlated with Instar (maximum number of instars within each species) (Tables 2 and 5). This correlation with Instar is likely a case where having the full ontogeny incorporated into an analysis would provide clarity on this odd result. In general, web builders tended to have less number of instars to adulthood than non-web builders or wandering spiders; exceptions in this study are the ambush thomisid and sit & pursue (Miller, Ament & Schmitz, 2014) philodromid species: Xysticus cristatus (C. Clerck, 1757) and Tibellus oblongus (C.A. Walckenaer, 1802), respectively (Table S1), which had relatively few instars. Some second instars possess the two MAP gland spigots observed in all species later in ontogeny. However, both X. cristatus and T. oblongus had no MAP gland spigots in the second instar and Metellina segmentata (C. Clerck, 1757) possessed only 1 MAP gland spigot (Table S1; Wąsowska, 1977; Yu & Coddington, 1990; Townley & Tillinghast, 2009). These three species had some of the lower maximum numbers of instars per species compared to the lycosoids, cribellate spiders and even viscous orb weavers (Table 2). Second instars across species possessed either none, one or two mAP gland spigots and this varied across foraging strategies and lineages. However, those species that possess both mAP gland spigots at the second instar were consistently the species with a higher maximum number of instars. This explains the significance of Instar as a predictor (Table S1).

Ancestral character estimation

We also conducted an ancestral character estimation on specific spigots on the posterior lateral spinneret whose potential homology have long been debated (Fig. 3; Peters, 1984; Peters, 1995; Griswold et al., 2005; Eberhard, 2010; Eberhard & Barrantes, 2015; R Carlson & CE Griswold, 1996, unpublished data; RE Alfaro, 2017, unpublished data). We added an additional five species, including two species sister to the Araneoids, and three cribellate species ancestral to both the Araneoidea and RTA clades (Fig. 3). As we previously described we constrained our rate matrix (Table 6) that we used to model substitution rates across branch lengths based on prior knowledge about historical possession of the flagelliform, pseudoflagelliform and modified spigots in each taxon and considering Dollo’s Law (Table 6) (Dollo, 1893; Platnick & Griswold, 1991; Griswold et al., 2005; Dimitrov et al., 2016; Wheeler et al., 2016).

The ancestral root of our phylogram of 27 species was more likely to have borne a modified spigot or none at all (likelihood: 0.365, 0.635, respectively) (Fig. 3, Results S5). Because this is not a full determination of the likelihood of a modified spigot or no spigot we cannot definitively determine the character state of the orb weaving ancestor at the node where Nicodamoidea + Araneoidea and the RTA clade split off (Fig. 3, Modified spigot likelihood: 0.398, None likelihood 0.602). It was unlikely to have possessed a pseudoflagelliform gland spigot (Fig. 3, Pseudoflagelliform likelihood: 0.000). We can, therefore, hypothesize that pseudoflagelliform and modified spigots are homologous structures and that modified spigots in the RTA clade likely are retained structures like those found in the sister and ancestral cribellate clades (Fig. 3). We cannot rule out that flagelliform gland spigots arose independently from modified spigots and thus cannot infer homology between this spigot and the others. We do not know the functionality of the modified spigots observed in D. tenebrosus and H. carolinensis. This would be useful to explore in the future, as the functionality of modified spigots in cribellate members of the RTA clade is the same as the pseudoflagelliform gland spigot (RE Alfaro et al., 2017, unpublished data). It is possible we are observing an intermediate stage of the loss of the modified spigot in D. tenebrosus and H. carolinensis.

Our results may change if we can incorporate ontogeny into the ACE analysis. For example, mimetids, which were ranked as having no spigots, do possess vestigial PLS spigots in the early instars (Table S1, Townley & Tillinghast, 2009). It is also possible that we are observing in real time the loss of the flagelliform gland spigot in this araneoid lineage. T. perfuga possess primordial modified spigot and flankers in second instars prior to them constructing webs in the third instar where functional spigots exist where the pre-spigots had been (Table S1; RE Alfaro, 2017, unpublished data). In most species, the final molt of the adult male leads to loss or nubbins (non-functional spigots) in all three: Modified, Pseudoflagelliform, Flagelliform (Table S1; Wąsowska, 1977; Yu & Coddington, 1990; Hajer, 1991; Townley & Tillinghast, 2009; Dolejš et al., 2014; R Carlson & CE Griswold, 1996, unpublished data). This coincides with the male abandonment of the web for an alternative lifestyle of wandering in order to find females.

Future work, as the methods of phylogenetic inference grow and progress, may yield different results if we incorporate this analysis with the triad of spigots, i.e., the triplet of MS plus flankers, PF and flankers and/or FL (flagelliform) plus AG (aggregate), or spigot associations, i.e., PC with the PF in cribellate orb weavers. In araneoids, the flagelliform gland spigot is flanked by two aggregate gland spigots and in T. perfuga the modified spigot is flanked by two spigots of unknown gland association. In Tengella, these flankers resemble AC gland spigots but in some other cribellate spiders, e.g., Matachia or Badumna (Griswold et al., 2005; figs 87 A, D) the flankers resemble paracribellar (PC) gland spigots. The color-coded ACE phylogram shows a trend of loss (black color) for all spigot types in approximately half of the tip species (Fig. 3). This coincides with a shift in foraging strategies, i.e., from webs to running, observed in these lineages (Table 2) and is consistent with the current hypothesis that an adaptive tradeoff between silk production and fecundity is driving spider evolution to foraging strategies that do not involve silk or web building (Blackledge, Coddington & Agnarsson, 2009). Given the recent conclusions of phylogenomics studies indicating a much more ancient orb weaving ancestor and the new sister relationship of cribellate orb weavers to the RTA clade, our ACE results do indicate that deeper and broader sampling across the spider tree of life is necessary (Bond et al., 2014; Fernández, Hormiga & Giribet, 2014; Garrison et al., 2016).

Increased taxon sampling would also improve the ancestral character estimation analyses performed for the independent variables of the PGLS analyses. Not surprisingly, given the results and discussion of the above analysis, the ancestral state at the root node of the phylogram of 27 species was Web builder for Strategy (foraging strategy; likelihood = 1.0), and specifically either funnel (likelihood = 0.812), sheet (likelihood = 0.084) or orb webs (likelihood = 0.105) (Figs. 9 and 10). Also, given the ancestral state for the singular spigot on the PLS as Modified spigot (MS), one may reason that the ancestor was also a cribellate silk user. We found this to be the case after performing ACE on the independent variable Silk (main type of silk expressed), and found the ancestral state at the root of the phylogram was 100% likely to be cribellate (Fig. 11). Finally, Type or variety of spigot types, was less conclusive and we found likelihood values for all states at the root node (Table S2). However, the most likely ancestral states were those possessing a cribellum, paracribellum (likelihood = 0.398) and possibly a pseudoflagelliform gland spigot (likelihood = 0.433) (Fig. 12). This may not actually conflict with our results for the ACE of the singular PLS spigot, i.e., that MS is ancestral, in that we determined that the MS and PF spigots are homologous structures. But we also found that the most likely state was no modified spigot (likelihood = 0.635), which would suggest the plesiomorphic state of cribellum + paracribellum, without any modifications of spigots on the PLS is the ancestral state at the root node of the ultrametric phylogram.

As the techniques for more complex phylogenetic comparative analyses improve, such as to allow for a time-series dataset with multiple values per species, we suspect that incorporating the entire picture of spigot ontogeny will lead to some interesting and novel inferences about silk evolution. By not incorporating the entire ontogeny, but “snapshots” of the adult female and second instars, important observations are missed by the analyses, such as loss and regain of AC gland spigots on the PLS in H. carolinensis, or presence of FL-AG triad spigots in early Mimetidae instars. This approach to understanding spigot ontogeny from a phylogenetic comparative perspective is novel and we can only build on our efforts from this study by growing the dataset to include deeper taxon sampling and working towards the capability of phylogenetic statistical analyses that can function to accommodate ontogeny datasets as whole units for each species.