Littorally adaptive? Testing the link between habitat, morphology, and reproduction in the intertidal sculpin subfamily Oligocottinae (Pisces: Cottoidea)

Phylogenetic framework

We constructed a phylogenetic framework using previously reported DNA sequences from all oligocottine species (sample size per species: 1–9 individuals, median 5) and several outgroups from the cottoid families (sensu Smith & Busby, 2014): Agonidae (n = 6 spp.), Cottidae (n = 1 sp.), Hexagrammidae (n = 1 spp.), Psychrolutidae (n = 11 spp.), and Rhamphocottidae (n = 1 sp.). These outgroup taxa are consistent with the most recent phylogenetic hypotheses of broader cottoid relationships (Knope, 2013; Smith & Busby, 2014). This dataset is accessible on Genbank (accession numbers KP826911–KP827632, see Table S1) and contains sequence data from eight molecular loci: one mitochondrial protein-coding locus (Cytochrome c oxidase, COI), two nuclear introns (exon-primed intron crossing (EPIC) locus 1777E10 and EPIC locus 4174E20) and five protein-coding nuclear loci (early growth response protein 1 (EGR1); mixed-lineage leukemia (MLL); patched domain-containing protein 1 (ptchd1); Rhodopsin; and Sushi, von Willebrand factor type A, and pentraxin domain-containing 1 (SVEP)). Multiple sequence alignments (MSAs) for each locus were generated in ClustalW (Larkin et al., 2007). Alignments were visually inspected, trimmed, and concatenated in Mesquite v3.2 (Maddison & Maddison, 2016). The best fitting model of molecular evolution for each locus was identified using the Akaike information criterion (Akaike, 1973; Posada & Buckley, 2004), with the model comparison routines implemented in MrModeltest v2 (Nylander, 2004).

The molecular dataset contains multiple representatives for each species, so we estimated a species tree using the multispecies coalescent model (Heled & Drummond, 2010) in BEAST v1.8.2 (Drummond et al., 2012). A species set was defined, based on the results of recent phylogenetic hypotheses (Knope, 2013; Smith & Busby, 2014), for the superfamily Cottoidea, which contains all taxa included in the dataset except Hexagrammos lagocephalus, which was treated as an outgroup. The species of each sampled individual in the dataset was assigned as a discrete trait. These species assignments were tested and validated in Buser & López (2015). For each locus, the model of molecular evolution yielding the lowest AIC value (as calculated in MrModeltest) was applied. The rate of molecular evolution was modeled as an uncorrelated lognormal relaxed clock (Drummond et al., 2006) and was unlinked across all loci. All tree models share a birth-death speciation tree prior with a piecewise linear and constant root population size model and a UPGMA starting tree. Four independent analyses were run for 200 million generations each and were sampled every 20,000 generations. MCMC logs were visualized using Tracer v1.6 (Rambaut et al., 2014) to determine convergence and an appropriate number of generations to discard as burn-in. Burn-in was removed and trees combined using LogCombiner v1.7.5 (Drummond & Rambaut, 2007). The phylogeny was pruned in the R statistical environment (R Core Team, 2015) using the “extract.clade” function from the package “ape” (Paradis, Claude & Strimmer, 2004) to only include members of the subfamily Oligocottinae plus the outgroup taxon Chitonotus pugetensis. An R script, “LitorallyAdaptiveScript.R”, detailing these commands and all other operations performed in R, along with all pertinent data (phylogeny, data matrix, etc.) is available in Supplemental Information 1.

Character coding

To understand the relationship between a species’ depth preferences and its size, reproductive habits, or scale patterns, we coded the following characters from previous studies and, where possible, verified our findings by examining museum specimens and/or collection data (summarized in Table 1):

1. Depth range (Continuous). Collection data for all specimens of each species of Oligocottinae and the outgroup taxon C. pugetensis were collated from museum records from the following natural history collections: University of Alaska Museum (UAM), University of British Columbia Beaty Biodiversity Museum (UBCBBM), University of Washington Burke Museum Fish Collection (UW), Oregon State University Fish Collection (OS), California Academy of Sciences (CAS), Natural History Museum of Los Angeles County (LACM), University of Michigan Museum of Zoology (UMMZ) and Scripps Institute of Oceanography Marine Vertebrates Collections (SIO). These records were accessed through institution-specific (UW, UBCBBM, CAS) or the multi-institutional database interfaces (all others) VertNet.org, Arctos.Database.Museum, and FishNet2.org (see Table S2 for all museum records analyzed). For each species, we extracted collection depth data from all museum holdings of adult specimens for which it had been recorded. Some collection depths are recorded as a range, in these cases, we used the maximum depth in the range. Where the collection depth and/or locality is described as “tide pool,” “intertidal,” etc., we assigned a collection depth of 0 m. To lessen the effects of outliers, we selected a depth range (i.e., minimum depth and maximum depth) for each species that includes 95% of museum collection depths (illustrated in Fig. 1). For the purposes of this study, we will refer to this depth range as the range where each species is “commonly” collected. To verify these depth ranges, maximum and minimum depth records for each species were cataloged and cross-examined from multiple sources Bolin, 1944; Miller & Lea, 1972; Eschmeyer, Herald & Hammann, 1983; Mecklenburg, Mecklenburg & Thorsteinson, 2002; see Table S3). Where these previously published depth maxima and minima disagree, we chose the median value for each. Many of these ranges include only imprecise descriptions such as “tide pools” and “intertidal areas”. In these cases, we assigned a minimum depth value of 0 m and a maximum depth value of 2 m.

2. Tide pool occupancy (Presence, absence). We noted which taxa were explicitly collected from tide pools in museum collection data, in previously published depth ranges, and in primary literature.

3. Length (Continuous). Maximum recorded length of each species was cataloged and cross-examined from multiple sources (Bolin, 1944; Miller & Lea, 1972; Eschmeyer, Herald & Hammann, 1983; Mecklenburg, Mecklenburg & Thorsteinson, 2002; Knope & Scales, 2013). Where sources disagreed, we used the median value.

4. Squamation (Presence, absence). For the purposes of this study, squamation is defined as any dermal ossification outside of the lateralis system. This includes scales, prickles, and scutes. The evolution of scale types in sculpins is poorly understood, but what is known suggests that the modified scales found in Oligocottinae may each represent an equal number of evolutionary steps away from the ancestral ctenoid scale type, with the latter not represented in any extant cottoid (Jackson, 2003). We therefore feel that in the context of this study it is unjustifiable to discriminate between scale types in oligocottines until further study indicates otherwise. Presence of squamation was coded from descriptions in the literature (Bolin, 1944; Begle, 1989; Mecklenburg, Mecklenburg & Thorsteinson, 2002; Jackson, 2003).

5. Enlarged genital papilla (Presence, absence). This character was coded directly from descriptions in the literature (Bolin, 1944; Mecklenburg, Mecklenburg & Thorsteinson, 2002).

6. Spermatozoon morphology (Oval, intermediate, slender). Character states were adapted from descriptions in the literature (Hann, 1930; Petersen et al., 2005; Koya et al., 2011). Slender sperm morphology is associated with internal insemination in many groups of fishes (Mattei, 1991). Petersen et al. (2005) confirmed this observation in Oligocottinae by demonstrating that spermatozoa with a slender-type morphology are active only in seawater that has been diluted to approximate the osmolality of ovarian fluid in these sculpins, while spermatozoa with oval-type morphology are active in both dilute and full-strength seawater. This suggests that slender-type spermatozoon morphology is indicative of obligate insemination, but oval-type morphology indicates the capacity for external mixing of gametes (i.e., spawning).

7. Copulation (Presence, absence). For the purposes of this study, copulation is defined as the transfer of sperm from a male into the ovary of a female. The presence of copulation, where known, was determined from descriptions found in the literature (Bolin, 1941; Morris, 1952; Morris, 1956; Hubbs, 1966; Misitano, 1980; Petersen et al., 2005; Abe & Munehara, 2009).

8. Parental care (Presence, absence). For the purposes of this study, egg guarding by one or both parents is considered parental care. The presence of parental care, where known, was determined from a review of behavioral descriptions from previous literature (Morris, 1952; Hubbs, 1966; Petersen et al., 2005; Abe & Munehara, 2009).

Character mapping and ancestral state reconstruction

To visualize the inferred evolutionary history of characters, we performed ancestral state reconstruction (ASR) of each character examined in this study. For discrete characters, we used maximum likelihood (ML) with a Markov k-state 1 parameter (Mk1) model of evolution (Lewis, 2001), implemented in Mesquite. The evolutionary history of continuous characters was inferred using ML in the R statistical environment with functions from the package “phytools” (Revell, 2012; see “Character mapping and ancestral state reconstruction” section in “LitorallyAdaptiveScript.R” in Supplemental Information 1).

Body shape analysis

Qualitative assessment of body shape in intertidal fishes has not only shown differences in the shape of some intertidal species compared to their subtidal relatives (e.g., Hypsoblennius spp.; Thomson & Lehner, 1976), but also that many groups of intertidal fishes (including intertidal sculpins) have converged on a small number of stereotypical body shapes (reviewed in Horn, 1999). While thought provoking, these observations have yet to be tested using quantitative methods. We used landmark-based geometric morphometrics to describe and compare the body shape of each species in this study and test for correlation between body shape and the depth at which each species occurs. Body shape data were collected from digital photographs of the lateral view of museum specimens of each species. Photography followed the phototank method of Sabaj Pérez (2009). To minimize the likelihood of introducing variation due to photographic artifacts (i.e., image distortion), lighting, distance to the subject, focal length, camera angle, and camera settings (e.g., aperture) were kept constant. We photographed 115 specimens, representing all 16 species in Oligocottinae plus the outgroup taxon C. pugetensis. Sample size per species ranges from 2 to 15 individuals, median six (Table 2). To capture overall head and body shape, landmarks were adapted from those described in previous studies of sculpin body shape (Strauss & Bookstein, 1982; Strauss & Fuiman, 1985). Preliminary analysis revealed a high frequency of distended stomachs and upturned caudal peduncles, presumably from prior feeding and preservation (respectively), so landmarks that appeared to be influenced by these variables were not included. Thirteen landmarks were ultimately used in this study (Fig. 2). We used tps-Dig2.2 (Rohlf, 2007) to locate the landmarks on each specimen from the digital photographs. To compare body shape across the group, landmark configurations were Procrustes superimposed using MorphoJ v1.06 (Klingenberg, 2011). The aligned landmark coordinates were used to calculate a covariance matrix on which we performed a principal component analysis (PCA) in MorphoJ. The number of significant principal component axes was calculated using the broken stick method (Frontier, 1976; Jackson, 1993; Legendre & Legendre, 2012), implemented with the “screeplot.cca” function in the R package “vegan” (Oksanen et al., 2017). The significant principal component axes were used to interpret overall shape variation and visualize the distribution of species in body shape morphospace. To visually check for evidence of morphological convergence or divergence, we projected phylogenetic relatedness into the principal component morphospaces and inferred states of each significant PC axis for each ancestral node (i.e., phylomorphospace analysis; Sidlauskas, 2008) using the “phylomorphospace” function in the R package “phytools” (Revell, 2012; see “Body shape analysis” section in “LitorallyAdaptiveScript.R” in Supplemental Information 1).

Depth correlation analysis

For both museum records and previously published depth ranges, preliminary results indicated that, while there is considerable variability in the maximum collection depth of each species in Oligocottinae, all species share a minimum recorded depth of zero meters. Given this invariability in minimum depth, we chose to use only maximum depth as our depth variable for regression analysis. We used phylogenetic generalized least squares regression (PGLS) implemented using the “gls” function in the R package “nlme” (Pinheiro et al., 2015) to test for a linear correlation between depth and each maximum length. To account for potential variability in trait evolution (e.g., Brownian motion, selection, etc.), we tested three alternate single-parameter correlation structures supplied in the R package “ape” (Paradis, Claude & Strimmer, 2004) in each of our regression models: a Brownian motion model with correlation due to phylogenetic relatedness represented by Pagel’s lambda (Pagel, 1994; Pagel, 1999), which we estimated using ML; a Brownian motion model with the rate of evolution (accelerated or decelerated) estimated using ML; and a single optimum (i.e., stabilizing selection) Ornstein–Uhlenbeck (OU) model (Felsenstein, 1988; Hansen, 1997) with the strength of attraction towards the optimum represented by alpha and estimated using ML. The best fitting model for each regression was determined by comparing AIC values. We also tested for a statistically significant difference in the average minimum or maximum depth between each state of the binary characters: presence of scales and presence of a genital papilla using phylogenetic analysis of variance (ANOVA) implemented using the “aov.phylo” function from the R package “geiger” (Harmon et al., 2007).

One of our hypotheses is that shallow-dwelling species show convergent morphology differing from that of their subtidal sister taxa. To test for convergent or divergent evolution of body shape, we used a stepwise model-fitting approach, “surface,” that detects shifts and convergence in phenotypic optima (Ingram & Mahler, 2013). In this approach, each optimum contributes a parameter to an OU process of evolution. The “surface” method finds the maximum-likelihood estimate of the number and location of phenotypic optima under the OU model and collapses similar phenotypic optima together if it improves the AIC score (Ingram & Mahler, 2013). Phenotypic convergence is indicated when independent lineages share a common optimum. These analyses were conducted in the R environment using functions from the package “surface” (Ingram & Mahler, 2013).

To visualize only the aspects of body shape that covary with depth, a partial least squares analysis (PLS) was conducted on a matrix of depth variables (minimum depth and maximum depth) and Procrustes-aligned shape variables (Rohlf & Corti, 2000). This analysis was conducted in MorphoJ and in R using functions from the package “geomorph v2.0” (Adams & Otárola-Castillo, 2013). In both cases, the significance of the covariance was tested using a permutation test with 10,000 iterations (see “Correlation with depth” section in “LitorallyAdaptiveScript.R” in Supplemental Information 1).

Character mapping and ancestral state reconstruction

The trimmed, concatenated MSA dataset spans 4,695 aligned nucleotide sites, containing 1,037 variable sites. The topology of the Bayesian maximum clade credibility (MCC) phylogeny produced herein is identical to the topology of the MCC phylogeny reported in Buser & López (2015), with similar levels of support for each clade (Fig. 3). As noted in Buser & López (2015), this topology is similar to that of other molecular-based phylogenetic inferences of Oligocottinae (i.e., Ramon & Knope, 2008; Knope, 2013), but has substantially higher support values (i.e., Bayesian posterior probability) for many of the inferred relationships. We will use the classification and taxonomy suggested by Buser & López (2015) for discussion of the interrelationships of oligocottine sculpins.

The outgroup taxon, C. pugetensis, rarely (if ever) occurs in intertidal areas (Fig. 1, Table 1, Tables S3 & S4). However, apart from L. hirundo, all the constituent species of Oligocottinae are regularly found in intertidal habitats and both museum records and published depth ranges include tide pools in the common collection depth or depth range data for all oligocottine species but L. hirundo (Fig. 1, Table 1, Tables S3 & S4). There is also explicit discussion of tide pool and intertidal occupancy for all oligocottine species except L. hirundo in the primary literature (Table S3). However, while the occupation of intertidal and subtidal habitats is often portrayed as an either/or scenario, there is considerable variation in the maximum depth at which each species occurs (Fig. 1, Table 1). Generally though, all oligocottine species occur at relatively shallow depths: none is commonly collected at depths greater than 55 m, most (12/16 spp.) are not commonly collected below 25 m (though there is some discrepancy between the museum collection data and the published depth ranges for A. corallinus and A. fenestralis), and four (published ranges) to seven (museum depth data) species are common only in very shallow (i.e., 2 m depth or less) habitats (Table 1). There is considerable disagreement between the museum collection data and the published depth range for A. notospilotus, C. acuticeps, C. analis, and L. hirundo. In each case, published depth ranges indicate a maximum depth that is >10 m deeper than the depths where these species have been commonly collected in museum holdings Tables S3 & S4). However, the depth ranges are otherwise largely congruent, and the differences between the two datasets are minimal (compare depth range values in Table 1 and Table S3). All remaining analyses show identical outcomes when using either the common museum collection depth data or the previously published depth range data for each species. Given the congruence of the datasets, the indistinguishable outcome of using one over the other, and the more verifiable nature of the museum collection records, we present the results of the remaining analyses using only the common museum collection depth range of each species.

Predictably, the ASR of minimum depth shows that the most recent common ancestor (MRCA) of Oligocottinae likely occurred in shallow habitats (ML estimate: 1 m; 95% confidence interval: 0 m, 2 m). Ancestral state reconstruction of tide pool occupancy shows that with extremely high proportional likelihood (0.9988) the MRCA of Oligocottinae occurred in tide pools. In fact, even the MRCA of the Leiocottus lineage was likely (0.9215 proportional likelihood) capable of living in tide pools (Fig. S1). Thus, the absence of tide pool occupation in L. hirundo likely represents a derived state. The ASR of maximum depth suggests that the MRCA of Oligocottinae occurred down to only moderate depths (ML estimate: 23 m; 95% confidence interval: 2 m, 44 m; see Fig. 4) and suggests that the habitation of only very shallow-water habitats (maximum depth = 2 m or less) seen in Oligocottus maculosus, O. rimensis, and O. snyderi and in all members of the subgenus Clinocottus (Blennicottus) represents a derived state (see Table 1, Fig. 4). However, given the uncertainty of the ML estimates of maximum depth at each node (Fig. 4), and the uncertain phylogenetic relationships of Blennicottus, Leiocottus, and Oligocottus lineages (Fig. 1), it is not possible to claim with confidence the number of transitions that may have occurred within the subfamily.

Maximum length shows no obvious relationship with depth in Oligocottinae (illustrated in Fig. 4). All species (including the outgroup) are relatively small (none longer than 250 mm), most (12/16 spp.) do not grow longer than 150 mm, and seven species do not grow longer than 100 mm (Table 1). The ASR of maximum length suggests that the MRCA of Oligocottinae was small (132 mm, 95% confidence interval: 83 mm, 182 mm), but like maximum depth, the uncertainty of the ML estimates at each node precludes additional inference into the diversification of this trait (Fig. S2).

Squamation is common among members of Oligocottinae but is completely absent in two clades: the subgenus Clinocottus (Blennicottus), and the clade in Oligocottus made up of O. maculosus, O. snyderi, and O. rubellio (Table 1, Fig. S3). Predictably, the ASR shows that the presence of squamation is the most likely state for the MRCA of Oligocottinae (proportional likelihood: 0.81), and that the absence of scales represents an independent loss of the trait in the MRCA of each of the above clades (Fig. S3). It should be noted that we inferred the evolution of squamation using an Mk1 model, which assumes that all potential changes in state are equally probable (Lewis, 2001). Given the lack of rigorous study of scale evolution in sculpins, specifying a more complex model is not warranted. However, it is our opinion that re-acquisition of squamation is an extremely unlikely evolutionary scenario in oligocottine sculpins (i.e., less probable than the loss of squamation) and thus the proportional likelihood of the presence of scales for the ancestral nodes in Oligocottinae should be taken as a conservative estimate.

An enlarged genital papilla is found in all but five species of Oligocottine sculpins (Fig. 5). The ASR shows with high proportional likelihood (>0.98) that this character was present at all ancestral nodes except those within the clade composed of the members of the genus Artedius. Within Artedius, an enlarged genital papilla was likely lost in the MRCA of the clade containing A. corallinus, A. fenestralis, A. lateralis, and A. notospilotus (Fig. 5). An independent loss of the enlarged genital papilla occurred in the oligocottinin species Orthonopias triacis.

The distribution of spermatozoon morphology follows a pattern similar to that of the enlarged genital papilla. Outside of the clade composed of members of the genus Artedius, slender-type spermatozoa are present in all Oligocottine sculpins and, with high proportional likelihood (>0.90), this is the state at all ancestral nodes (Fig. 5). Within Artedius, an oval-type spermatozoon likely evolved in the MRCA of the clade containing A. corallinus, A. fenestralis, A. lateralis, and A. notospilotus. Artedius harringtoni possess an intermediate spermatozoon morphology and is the only oligocottine to do so (Fig. 5). With one notable exception, this shows that all species with a slender-type spermatozoon morphology (which in other species is known to function only in ovarian fluid) also possess an enlarged genital papilla, which is presumably used in copulation. The exception to this observation is the species Orthonopias triacis, which does not possess an enlarged genital papilla, but does possess a slender-type sperm morphology.

Though not known for all species in Oligocottinae, the distribution of copulatory behavior closely follows that of spermatozoon morphology, where copulating species possess either slender or intermediate-type spermatozoon morphology, and non-copulating species possess only oval-type spermatozoon morphology (illustrated in Fig. 5; see also Fig. S4). The ASR of this character shows that, with high proportional likelihood (>0.95), presence of copulation is the likely state for all ancestral nodes outside of the clade composed of the members of the genus Artedius. Within Artedius, copulation was likely lost in the MRCA of the clade containing A. corallinus, A. fenestralis, A. lateralis, and A. notospilotus. This finding reinforces the observations of Mattei (1991) and Petersen et al. (2005) who each show that slender-type sperm morphology is indicative of copulating species. Once again, Orthonopias triacis presents a noteworthy case as there is evidence of copulation for the species and the species possesses a slender-type spermatozoon morphology, yet the species lacks an enlarged genital papilla or other known intromittent organ (Fig. 5).

The presence or absence of parental care has been described in less than half of all oligocottine species (6/16 spp.), but follows a similar pattern to those seen in other reproductive characters in the group, where members of the genus Artedius tend to differ from all other species. In this case, parental care is observed only in members of Artedius (Fig. 6). The ASR shows that parental care was likely present in the MRCA of Artedius, while a lack of parental care is the most likely state for the MRCA of the tribe Oligocottini. However, given the substantial amount of missing data for this trait, the ASR is subject to change with the addition of new observations.

Body shape analysis

Observed body shape variation was captured by two significant principal components, which cumulatively describe 70% of the total variance. We used an outline of a specimen of Clinocottus analis to visualize shape change represented by each PC axis in MorphoJ (Fig. 7). Principal component (PC) 1 (52% of total observed variance) describes antero-posterior elongation/compression of the head and mouth as well as the relative size of the eye. Principal component 2 (19% of total observed variance) captures dorso-ventral elongation/compression of the body, the shape of the eye, and the slope of the snout.

There is no clear evidence of a consistent relationship between the minimum depth or maximum depth of a species and its morphology, nor evidence of morphological convergence among shallow or deep-ranging species (Fig. 8). Species with deeper ranges appear to be constrained to a common morphospace, while species that inhabit only shallow depth ranges (e.g., Oligocottus spp., Clinocottus (Blennicottus) spp.) appear to occupy novel and distinct areas of morphospace (Fig. 8). This observation is supported by the results of the “surface” analysis, which inferred three optima for body shape in the morphospace described by the significant PC axes: one for members of the genus Oligocottus plus Orthonopias triacis, one for members of the subgenus Clinocottus (Blennicottus), and a third that is occupied by all other oligocottines (Fig. S5). However, there is no clear pattern in terms of the direction of the divergence in morphospace of these taxa. Interestingly, two of the optima are occupied almost entirely by taxa that are found exclusively in shallow water (i.e., Clinocottus (Blennicottus) and all but one species of Oligocottus), while the remaining optimum is made up almost exclusively of deeper-ranging taxa (the exception being C. (O.) acuticeps. Here again Orthonopias triacis is remarkable in that it is a deeper-ranging species, but appears to be drawn to the phenotypic optimum occupied otherwise exclusively by the genus Oligocottus.

Depth correlation analysis

No morphological, reproductive, and body shape variables examined in this study show a significant correlation with maximum depth. The phylogenetic generalized least squares regression showed no significant linear correlation between depth and maximum length, and the phylogenetic ANOVA showed no statistically significant difference in average maximum depth between the binary character states of: presence of squamation (p-value > 0.15) or presence of an enlarged genital papilla (p-value > 0.70). Likewise, the results of the PLS analysis were identical in MorphoJ and in R, and failed to show a statistically significant correlation between body shape and depth range (p-value > 0.05). However, the analysis did reveal a few interesting observations: species that are restricted entirely to shallow water tend to display greater morphological diversity than those that range into deeper water, but generally have smaller eyes, smaller mouths, terminal placement of the mouth, and more robust bodies, while deeper-ranging species tend to have larger eyes, larger mouths, subterminal placement of the mouth, and dorso-ventrally compressed bodies (Fig. 9).

Squamation, length, and depth

The subfamily Oligocottinae should be thought of as a clade of intertidal-occurring fishes and the ability to live in intertidal depths and specialized intertidal habitats such as tide pools is likely the ancestral state of the group. However, this ability to live in tide pools does not preclude residency in other habitat types within the same species, as many of the extant and ancestral species are capable of living in a variety of depths in addition to the intertidal ones. Thus, the diversification of Oligocottinae should be characterized asoccurring within a habitat range that includes both subtidal and intertidal habitats. This may explain the general lack of correlation between depth the other characters examined in this study.

Small maximum size and a reduction in scales have been noted as common features of intertidal fishes by previous authors (Gibson, 1982; Knope & Scales, 2013), and while we found no evidence to support these hypotheses within Oligocottinae, oligocottines as a whole may in fact offer support. All oligocottines are small (none longer than 250 mm), and all show a reduction in scales when compared to a “typical” scaled member of Cottales, such as Oxylebius pictus (see Jackson, 2003). Outside of the lateral line, the most heavily scaled oligocottines possess only modified scales in a narrow band on the body along each side of the dorsal fins and on the dorsal surfaces of the head and caudal peduncle (e.g., Orthonopias triacis; see descriptions in Bolin, 1944; Jackson, 2003). Others possess highly reduced scales in the form of prickles (e.g., Clinocottus (Clinocottus) analis; see description in Bolin, 1944). It is presumed that the primary reason scales are often reduced in intertidal fishes is to promote cutaneous respiration, which is dependent on well vascularized skin that is free from obstructions (Feder & Burggren, 1985; Martin & Bridges, 1999). Within Oligocottinae, it may simply be the case that the highly-reduced scales seen in the group do not cover enough surface area to interfere with cutaneous respiration in a meaningful way. Or perhaps that cutaneous respiration is restricted to only certain areas, such as the head (as seen in Coryphoblennius galerita, see Zander, 1972), or buccal chamber (reviewed in Bridges, 1993). This may explain the presence of scales within the group, and the fact that even the most heavily scaled members of Oligocottinae occur intertidally and in tide pools (Fig. 2). Interestingly, the one oligocottine species that does not occur intertidally (L. hirundo) possess only a few highly reduced scales in the form of a scattered patch of prickles located just posterior to the base of the pectoral fins (see Jackson, 2003). We interpret this as evidence that a reduction in scales in an ancestral condition for Leiocottus, and possibly Oligocottinae as a whole.

The evolutionary history of maximum size, depth range, and the presence of squamation all support the idea that the subfamily Oligocottinae is a primitively intertidal group. All species in this subfamily display conditions of these traits that are common in intertidal taxa, and all but one species are commonly found intertidally and explicitly in tide pools. These findings support those of previous studies which found Oligocottinae to represent an intertidal radiation (Ramon & Knope, 2008; Knope & Scales, 2013).

Body shape and depth

While the body shape of deeper-ranging species is relatively conserved, groups that occupy only intertidal areas appear to be exploring novel areas of morphospace (Fig. 8, Fig. S5). This is most pronounced in members Oligocottus and the subgenus Clinocottus (Blennicottus), where the two groups each occupy a morphological optimum that is distinct from each other and from all but one other oligocottine (i.e., the enigmatic O. triacis). The morphological distinctiveness of Clinocottus (Blennicottus) is readily apparent, even to the casual observer, with antero-posteriorly compressed heads and highly robust bodies and fins (see illustrations in Bolin, 1944). The constituent species of this subgenus occur most abundantly in steep, rocky habitats with high wave exposure and are often the only oligocottine species found in these areas (T. Buser pers. obs.; Green, 1971; Yoshiyama, 1981; Yoshiyama, Sassaman & Lea, 1986; Mgaya, 1992, see also descriptions in Eschmeyer, Herald & Hammann, 1983; Lamb & Edgell, 1986; Mecklenburg, Mecklenburg & Thorsteinson, 2002). The blunt heads and short, stocky bodies of these species are also seen in other intertidal fishes occupying similarly exposed, rocky habitats and may reflect a common evolutionary response to the physical demands of living in such areas (Kotrschal, 1988; Kotrschal, 1989; Thomson, Findley & Kerstitch, 2010, reviewed in Kotrschal, 1999). If the diverse morphologies seen among intertidal specialist clades are reflective of their respective habitat partitions, it may also be the case that the constrained morphologies seen in deeper-ranging taxa reflect a kind of stabilizing selection of generalist traits that are optimal for occupying a comparatively wide variety of habitats. The relatively high diversity of morphotypes seen in shallow vs deep-ranging species may mask morphological traits that are in fact associated with depth, as is suggested by the results of our PSL analysis (Fig. 9). However, the lack of statistical significance of this trend warrants caution on interpretation of this finding until the question can be revisited with additional taxa.

Reproduction

While reproductive characters do not show any correlation with depth range in Oligocottinae, the evolution of these traits in the subfamily may offer new insight into the evolution of reproductive modes in cottoids. Copulation in oligocottines is associated with an enlarged genital papilla and a slender-type spermatozoon morphology. These traits are broadly distributed in Oligocottinae and are likely the ancestral state of the subfamily (Fig. 5). Importantly, the absence of copulation and associated traits in most members of the genus Artedius represents a loss and is thus a derived state. This finding runs counter to previous hypotheses of the evolution of reproductive modes in sculpins, which interpreted the seemingly scattered distribution of copulation in cottoids as indicative of parallel or convergent evolution of copulation from non-copulating ancestors (Abe & Munehara, 2009; Muñoz, 2010). Under this paradigm, Petersen et al. (2005) suggests that the ability of the oval-type sperm morphology (uniquely capable of functioning well in seawater and ovarian fluid, seen in non-copulating members of Artedius) to function in ovarian fluid represents a derived condition and concludes that the presence of this trait in most members of Artedius represents an evolutionary step towards copulation in the group. We conclude the opposite of Petersen et al. (2005), and suggest that rather than the sperm’s ability to function in ovarian fluid, it is in fact the sperm’s ability to function in seawater that is a derived state and this, along with the loss of an intromittent organ, represents an evolutionary step away from copulation within Artedius.

A reduction or loss of the enlarged genital papilla is seen in other oligocottines as well. Critically, however, these species maintain a slender-type spermatozoon morphology and, where known, copulation. For example, while most oligocottines possess genital papillae that are quite large and robust, males in the genus Oligocottus possess papillae that are uniquely small, gracile, and thread-like (Fig. S6). Also unique to the males of this genus are modifications of the anterior portion of the anal fin (Fig. S7) which, where known, is used for grasping females during copulation (O. snyderi; Morris, 1956). It is possible that the added security and stability during copulation provided by the prehensile anal fin rays has rendered the large genital papilla seen in other oligocottines redundant. The other example of a reduction in the size of the male genital papilla is seen in Orthonopias triacis, where males lack an intromittent organ altogether, yet also possess slender-type sperm morphology and are known to copulate. Males of this species possess enlarged pelvic fins that face inwardly “palm to palm,” and project postero-ventrally from a “pedunculated” base (Bolin, 1944). Perhaps these highly modified, sexually dimorphic pelvic fins are used in a grasping manner that, like in Oligocottus, is used to aid in copulation and has eliminated the need of a large, robust male genital papilla. Copulation without the use of an intromittent organ is seen in at least one other member of Cottales, the sea raven (Hemitripterus villosus). In this species, males are not known to possess any putative grasping organs. Instead, the female everts her genital tract, which is covered in mucus, and the male ejaculates onto it, whereby the sperm become entrained in the mucus and enter the female when she inverts her genital tract (Munehara, 1996).

Our results show that while the presence of an enlarged male genital papilla is a likely indicator of copulation, the absence of an intromittent organ does not necessarily indicate the absence of copulation. Furthermore, our results show that non-copulating species may evolve from copulating ancestors. Given the widespread distribution of copulation and/or an enlarged genital papillae within Cottoidea (Abe & Munehara, 2009; Muñoz, 2010), we suggest that copulation and associated traits may have evolved much earlier in cottoids than has been previously estimated. Perhaps the seemingly disparate distribution of copulation in cottoids is not due to many independent evolutions of copulation, but rather to a single early evolution of copulation and multiple subsequent losses of the trait. Given the suite of complex physiological and behavioral traits associated with copulation in sculpins (e.g., internal gamete association with delayed fertilization, see Munehara, Takano & Koya, 1989; Munehara, Takano & Koya, 1991; Munehara et al., 1997; Petersen et al., 2005), the independent loss of copulation by certain lineages would, in our opinion, be a far simpler explanation for the modern distribution of the trait than would the independent evolution of copulation and all associated characters.

Parental care

Like other reproductive traits, the distribution of parental care in Oligocottinae does not appear to be related to the distribution of depth ranges. Rather, only members of the genus Artedius display parental care. Strong phylogenetic signal of parental care has been reported for other groups as well (reviewed in Coleman, 1999), but this does not provide a satisfying explanation for why Artedius differs from all other oligocottines in this trait. It is interesting to note that parental care shows an almost inverse distribution to oval-type sperm morphology (i.e., obligate copulation; see Figs. 5 and 6), but the limited sample size and degree of missing data for parental care make this a tenuous connection. Many other non-copulating sculpins also display parental care (e.g., Hemitripterus spp., Enophrys bison, Myoxocephalus spp., Cottus spp.), but this trait is also seen in some copulating species, including Artedius harringtoni (Abe & Munehara, 2009). Clearly, more research is needed to better understand the evolution of parental care in cottoids, and its relationship with other aspects of their complex reproductive biology.