The evolution of jaw protrusion mechanics is tightly coupled to bentho-pelagic divergence in damselfishes (Pomacentridae)

Transitioning to a different trophic niche typically requires or is accompanied by a change in the functional morphology of feeding (Christiansen and Wroe, 2007; Liem, 1980b; Westneat, 1994). Protrusile jaws have evolved multiple times among fishes (Staab et al., 2012; Wainwright et al., 2015). Two of the most successful vertebrate lineages, the Acanthomorpha (∼17,000 species) and Cypriniformes (∼3200 species), are composed of fish species that have rapidly transitioned between forms with highly protrusile upper jaws and those that exhibit little to no jaw protrusion (Hernandez and Staab, 2015; Hulsey et al., 2010; Staab et al., 2012; Wainwright et al., 2015; Bellwood et al., 2015). Both acanthomorphs and cypriniforms possess jaw arrangements that allow simple shape changes to either enhance or reduce jaw protrusion, with changes in the length of the ascending arm of the premaxilla having particular importance (Hernandez and Staab, 2015; McGee et al., 2015b; Rice et al., 2008; Staab et al., 2012). Both clades also occupy an extremely large number of feeding niches (Hernandez and Staab, 2015; Wainwright et al., 2015). If this diversity is largely a product of having protrusion mechanisms in which small morphological changes produce adaptive functional shifts, then these jaw mechanisms can be regarded as highly evolvable biomechanical systems (Pigliucci, 2008). This is particularly true if such changes are likely to arise through normal developmental variation in morphogenesis.

Aquatic environments have two primary sources for food: the water column (pelagic feeding) and the substrate (benthic feeding). Animals may obtain sustenance from either realm exclusively or they may occupy trophic niches that lay along a continuum between these extremes. One of the most common patterns of evolution among aquatic animals is transitioning between feeding niches that lie at different points along this bentho-pelagic niche axis, and such diversification has arisen repeatedly in molluscs, crustaceans, annelids, pinnipeds, elasmobranchs and bony fishes (Bracken et al., 2009; Cooper et al., 2010; Jones et al., 2013; Lindgren et al., 2012; Regier et al., 2010; Struck et al., 2015; Wilga et al., 2007).

We investigated whether the evolution of jaw protrusion ability has been linked to divergence along the bentho-pelagic axis in an adaptive radiation of acanthomorph fishes: the damselfishes (Pomacentridae). The Pomacentridae are a successful lineage of nearshore reef fishes (399 extant species; Eschmeyer and Fricke, 2016) that are one of the dominant vertebrate groups on coral reefs (Bellwood and Hughes, 2001; Allen, 1991; Bellwood et al., 2016). These fishes are highly amenable to this type of study because the evolution of their cranial morphology has tracked repeated transitions between benthic and pelagic feeding niches (Cooper and Westneat, 2009; Frédérich et al., 2013; Olivier et al., 2016). Their value as an experimental system is further enhanced by the relative ease with which they may be captured in the wild or purchased through the aquarium trade, the readiness with which most damselfishes feed in aquaria and the amount of published information on their diets. Although damselfish feeding mechanics have undergone rapid evolution, this has not resulted in the exploitation of a large number of food resources, but has instead produced a pattern of ‘back and forth’ shifts between only three primary trophic states: pelagic feeding on plankton, benthic feeding on algae (which may include significant feeding on detritus; Wilson and Bellwood, 1997) and bentho-pelagic omnivory, which involves feeding on a mixture of algae, plankton and small benthic invertebrates (Cooper and Westneat, 2009; Frédérich et al., 2013). This pattern of rapid ecomorphological evolution in conjunction with repeated invasions of the same niches, but without diversification into new niches, has been described as a reticulate adaptive radiation (Cooper and Westneat, 2009).

Quantitative studies of functional diversity are far less common than studies of morphological diversity, and this can be at least partially attributed to the extensive time required to characterize the functional abilities of large numbers of taxa (Wainwright, 2007). Because multiple morphologies can yield similar performance capabilities, inferring functional properties from morphological data alone is problematic (Wainwright, 2007; Wainwright et al., 2005) and combined studies of anatomical and functional diversification are necessary if we are to accurately map form to function relationships. Furthermore, the use of phylogenetic comparative methods to analyse functional data has been limited, and this impairs the validity of many of the statistical methods that have been used to examine functional diversification. Although it has been demonstrated that differences in jaw protrusion mechanics between closely related species are associated with differences in trophic ecology (Hernandez and Staab, 2015; Holzman et al., 2008a; McGee et al., 2015b; Rice et al., 2008), we know of only a single study that corrected for relatedness in order to identify which aspects of functional morphology have evolved in correlation with jaw protrusion (Hulsey et al., 2010). No previous studies have used phylogenetic comparative methods to: (1) determine whether fishes in different trophic guilds possess different protrusion abilities, or (2) characterize which evolutionary models best fit patterns of divergence in protrusion ability among the members of a lineage.

We performed combined analyses of skull form and function using 10 damselfish species that include multiple examples of convergence on all three primary pomacentrid trophic states, and analysed our data using phylogenetic comparative methods (Fig. 1, Table 1). We collected high-speed video of all specimens as they fed from the water column when jaw protrusion is most pronounced, as opposed to feeding via biting on attached food items when protrusion is less evident or absent. Feeding from the water column is at least occasionally employed by all damselfishes because even predominantly benthic-feeding species are known to do so during times of high pelagic food abundance (McCormick, 2003; Pratchett et al., 2001; Westneat and Resing, 1988).

We compared the trophic morphology of these 10 species to those of an additional 40 damselfishes from all 29 pomacentrid genera in order to: (1) confirm that the filmed species represent wide coverage of the anatomical diversity of the Pomacentridae and (2) better describe those skull morphologies that enhance protrusion ability. We then performed evolutionary analyses of form, function and diet data in order to test the following hypotheses: (1) damselfishes in trophic niches from different points on the bentho-pelagic spectrum have significant differences in jaw protrusion ability; (2) jaw protrusion ability has evolved in correlation with a suite of additional morphological and functional traits associated with feeding; (3) pomacentrids capable of extensive jaw protrusion exhibit higher levels of functional integration between upper and lower jaw movement than other species; and (4) patterns in the diversification of damselfish feeding mechanics are best described by evolutionary models of adaptation to three trophic niches: planktivory (pelagic feeding), omnivory (bentho-pelagic feeding) and herbivory (benthic feeding).

All aspects of this study, including fish euthanasia, were performed in adherence with Washington State University IACUC protocol 04285.

Specimens of Amphiprion frenatus Brevoort 1856, Amphiprion ocellaris Cuvier 1830, Chromis cyanea (Poey 1860), Chrysiptera cyanea (Quoy and Gaimard 1825), Dascyllus aruanus (Linnaeus 1758), Lepidozygus tapeinosoma (Bleeker 1856) and Pomachromis richardsoni (Snyder 1909) were obtained from the pet trade. Specimens of Chromis viridis (Cuvier 1830), Pomacentrus moluccensis Bleeker 1853 and Stegastes nigricans (Lacepède 1802) were collected using dip nets and barrier nets from reefs around the Lizard Island Research Station on the northern Great Barrier Reef, Australia. Preserved specimens from an additional 40 species were obtained from The Field Museum (Chicago, IL, USA) and used in morphological analyses (see below). Our taxonomy follows Cooper and Santini (2016).

After the feeding trials, the heads of all specimens were dissected to expose the functional morphology of the oral jaws (Fig. 2) (see Cooper and Westneat, 2009 for further details). Fishes were euthanized in accordance with approved IACUC protocols, formalin fixed until rigid, leached of formalin in tap water and then stepped over into 70% ethanol. Photographs of all dissected heads were taken in lateral view with the mouth closed using either a Nikon Coolpix S8200 digital camera or an Olympus DP25 digital camera interfaced with an Olympus SZ61 dissecting microscope. A scale bar was included in each photograph. One C. viridis specimen did not fix properly and was not photographed. In our discussions of jaw muscle morphology we follow the anatomical nomenclature used by Datovo and Vari (2013). The pars malaris, pars rictalis and pars stegalis divisions of the adductor mandibulae (the major biting muscle for most teleosts) are synonymous, respectively, with the A1, A2 and A3 nomenclature for these divisions established by Winterbottom (1973). Our references to cranial bone morphology follow Barel et al. (1976). After dissected whole heads had been photographed, we removed a premaxillary bone (upper jaw) and one side of the lower jaw (i.e. mandible, which consists of left and right articular and dentary bones in teleosts) from each specimen and photographed these skeletal elements in lateral view.

The coordinate locations of 18 anatomical landmarks (LM) of functional importance to feeding (Fig. 2, Table S1) were obtained from digital images of dissected heads using the program tpsDIG2 (http://life.bio.sunysb.edu/morph/). We used an outline-based semi-landmark approach to obtain shape data for the premaxillae and lower jaws. The programs tpsUtil and tpsrelW were then used to superimpose these semi-landmarks using a chord-distance (Procrustes distance) based ‘sliders’ method (http://life.bio.sunysb.edu/morph/). Overall head shape, premaxillary shape and lower jaw shape were analysed separately using principal components analyses (PCA).

For all shape analyses we used the program CoordGenMac7a (http://www3.canisius.edu/~sheets/imp7.htm) to calculate the Procrustes mean shape for each species. Once mean LM configurations were calculated, a second Procrustes transformation was performed using the mean shapes of all 10 species to remove differences in size or orientation from the coordinate data. We then performed a phylogenetic PCA (pPCA) on the transformed skull shape data in order to reduce the LM data to a smaller number of independent, orthogonal axes, and correct for different degrees of relatedness among species. We performed pPCA using the phyl.pca function in phytools (Revell, 2012).

A second head-shape pPCA analysis was performed using the filmed specimens plus LM data from an additional 40 species (Table S1) representing all damselfish genera (two to three specimens per species). This was done in order to determine how the head shape diversity of the 10 filmed species compared with the extant head shape diversity that has evolved within the damselfish lineage as a whole. For this larger dataset, only 16 head landmarks were analysed. Data for the pars stegalis (A3) division of the adductor mandibulae were omitted because the LMs for this muscle (LMs 7 and 15) were not visible on all images of the additional 40 species.

All fishes were acclimated to glass aquaria and filmed feeding on commercial food pellets. Fishes were filmed in lateral view at 250 frames s⁻¹ during feeding strikes using a Redlake high-speed video camera (Redlake MASD, San Diego, CA, USA) or an Edgertronic monochrome high-speed video camera (Sanstreak Corp., San Jose, CA, USA). Pellets were broken into pieces for smaller fishes and we attempted to maintain a consistent pellet-to-fish size ratio among species with the largest species (Chromis cyanea and Stegastes nigricans) receiving whole pellets. We selected three to four clear strike sequences for each of the three specimens per species.

On every other video frame, the coordinate locations of the food pellet and seven anatomical landmarks on each fish (Fig. 2) were recorded using the program tpsDIG2 (http://life.bio.sunysb.edu/morph/). The ImageJ software program (Schneider et al., 2012) was used to collect kinematic data from changes in LM positions. The following kinematic variables were measured (see Fig. 2 for identification of LM): distance moved by the fish (linear distance between points 2 and 8), distance moved by the food pellet (linear distance between points 7 and 8), distance between fish and pellet (linear distance between points 3 and 7), gape (distance between points 3 and 4), gape angle (angle 3, 5, 4), jaw protrusion (distance between points 2 and 3) and cranial elevation (angle 1, 2, 6). Variables likely to be affected by differences in the size of the fish, such as gape, jaw protrusion and velocity, were scaled to the fish's standard length. Velocity and acceleration were calculated as first and second derivatives of time and distance to the fixed point (LM 8). These data were used to calculate the kinematic variables listed in Table 2. Kinematic plots were used to track and compare changes in individual variables over time and all motion sequences were aligned to the time at which the food pellet first entered the mouth (time zero). In each feeding event, the distance (D) traveled by the fish and the distance moved by the food pellet (as a result of suction produced by the fish) were used to calculate the ram-suction index (RSI) following Norton and Brainerd (1993): RSI=(D_fish−D_pellet)/(D_fish+D_pellet). This measure specifically examines ‘body ram’ (movement of the body toward a food item) and not ‘jaw ram’ (movement of the jaw toward a food item; Liem, 1980a).

In acanthomorph fishes, upper jaw protrusion is driven by lower jaw abduction via a rotational maxillary bone in the upper jaw that transmits abduction of the lower jaw (mandible) into anterior motion of the premaxillae (Westneat, 1990). The functional integration of premaxillary protrusion with lower jaw abduction was calculated as the coefficient of determination (r² value) of a quadratic relationship fitted to these variables (Wainwright et al., 2008), with lower jaw abduction angle as the independent variable.

For phylogenetic comparative analyses we used a trimmed consensus tree (Fig. 1) (Cooper et al., 2009). We used the contMap function in the R package phytools (Revell, 2012) to produce contour map phylogenies that depict estimates of the evolution of all 18 variables listed in Table 2. A phylogenetic generalized least squares (PGLS) analysis was used to examine the relationships among our kinematic (N=10), shape (N=8) and integration data (N=1). Given the shared evolutionary histories of the damselfishes, more closely related species are assumed to exhibit more similar trait values (Revell, 2010). As a result, more similar residuals are produced from a least squares regression line. PGLS takes into account the expected covariance structure of these residuals, and generates new slope and intercept estimates to account for interspecific autocorrelation due to phylogeny (Symonds and Blomberg, 2014). To account for phylogenetic structure, we first extracted the expected covariance under a Brownian model from our tree using the corBrownian function in ape (Paradis et al., 2004). We then used the gls function in the nlme package to perform the PGLS analysis.

We used a multivariate model-fitting approach with the R package mvMORPH (https://cran.r-project.org/web/packages/mvMORPH/index.html) to compare the fit of multiple models of trophic evolution to the following kinematic and morphological variables: MJP, RSI, LDA, A1MA and A2MA (see Table 2 for the key to all variable abbreviations). Model fitting of large numbers of variables is problematic, so we chose this reduced set of five variables based on their strong association with the functional morphology of fish feeding. We examined the fit of three single rate models: Brownian motion (BM), early burst (EB) and a single-peak Ornstein–Uhlenbeck (OU 1) model. Support for the BM model would suggest that trophic functional morphology is uniformly increasing over time. Support for the EB model would suggest that most of the disparity present in the trophic functional morphology of modern damselfishes was partitioned early on in their evolutionary history (Harmon et al., 2010). Support for the OU 1 model would suggest there is a single, optimal combination of head morphology and bite kinematics for all of the Pomacentridae. The single-peak model was compared with five multi-peak OU models (a single two-peak model; three models with three peaks; and a single four-peak model) in which each of the 10 species were assigned to trophic categories based on published diet data (Table 1). To account for uncertain character histories in our five multi-peak models, we used the Stochastic Mutational Mapping on Phylogenies (SIMMAP) tool (Bollback, 2006) from phytools (Revell, 2012). We produced 500 simulated character history trees for each of the five multi-peak models and ran our multivariate trait data over each SIMMAP tree. We used the second-order Akaike’s information criterion (AICc), which corrects for small sample sizes, to select among the best evolutionary models. We calculated a mean and 95% confidence interval AICc score for each of the five multi-peak models. The best-fitting model is determined by the lowest AICc score and is favored over any other model if the difference in AICc score is greater than two units (Burnham and Anderson, 2002).

Performing evolutionary model analysis on a sample with few species can result in high error rates depending on the structure of the data (Boettiger et al., 2012; Goolsby, 2016). To determine whether we had the statistical power to detect an effect in our model comparison study, we performed a phylogenetic power analysis using the phylocurve package in R (Goolsby, 2016).

We used the R package geomorph (Adams and Otárola-Castillo, 2013) to perform phylogenetic ANOVA (pANOVA) tests in order to determine whether the kinematic, morphological and integration variables we measured differed among damselfishes in different trophic classes. We used the trophic class assignments that matched those for the best-supported evolutionary model (Tables 1, 3). The pANOVA test statistic was calculated from the data and compared with a null distribution generated via 1000 simulations of new dependent variables determined from a single rate matrix on the phylogenetic tree. For those variables for which the initial pANOVA test returned significant results, we performed post hoc phylogenetic pANOVAs that compared members of each of the trophic classes to one another: planktivores versus omnivores; planktivores versus herbivores; and omnivores versus herbivores.

To quantify the degree of convergence on specialized trophic ecomorphologies, we used the convnumsig function from the convevol R package (Burd et al., 2014; https://cran.r-project.org/web/packages/convevol/index.html). Convnumsig performs 1000 evolutionary simulations along our pomacentrid phylogeny using parameters derived from the observed cranial shape data based on the first four pPC scores. We defined taxa residing in derived pelagic and benthic regions of morphospace (see Fig. 3A) and calculated the number of convergent transitions into each of those regions.

The 10 species filmed during feeding represent wide sampling of the damselfish radiation in terms of both phylogenetic and morphological diversity (Figs 1, 3). The first principal component derived from the pPCA (PC1) of the 50 species dataset accounts for 44.27% of the total shape variation and describes differences in eye size, mouth size, orientation of the mouth (terminal to subterminal), height of the supraoccipital crest, size of the pars malaris (A1) and pars rictalis (A2) divisions of the adductor mandibulae, dorso-ventral width of the jaw bones, proximity of the insertion of the pars malaris to the maxillary–palatine joint and the length of the ascending arm of the premaxilla. The shape variation associated with PC1 is well represented by the shape differences depicted in Fig. 3B.

Within the 50 species morphospace, the head shapes of the filmed species Chromis cyanea, Chromis viridis and Lepidozygus tapeinosoma clustered within a distinct region of damselfish head-shape space (Fig. 3) occupied by fishes that feed primarily from the water column (Emery, 1983; Gerber and Marshall, 1974; Mohsin et al., 1986; Morris, 1984; Randall, 1967; Sano et al., 1984). Because the diets, head shapes and feeding kinematics (see below) of these species were largely distinct from those of other pomacentrids, we will refer to these species as ‘derived pelagic feeders’. The head shapes of an additional five species also fell within the derived pelagic feeder area of morphospace (Fig. 3A) and there was a significant level of convergence on head shapes that lie within this region (convnumsig; P<0.01). The results of evolutionary convergence analysis also identified five instances of convergence on derived pelagic-feeding head shapes (Figs 1, 3).

The head shape of Stegastes nigricans clustered with another group of damselfishes with distinct head shapes that feed primarily from the benthos (Grove and Lavenberg, 1997; Hobson, 1974; Letourneur et al., 1997; Montgomery, 1980; Randall, 1985; Sano et al., 1984; Sikkel, 1995; Wilson and Bellwood, 1997). Since the diets, head shapes and feeding kinematics (as represented by S. nigricans) of these species were largely distinct from those of other pomacentrids (see below) we will refer to these species as ‘derived benthic feeders’. The head shapes of an additional five species also fell within the derived pelagic feeder area of morphospace and there was a significant level of convergence on head shapes that lie within this region (convnumsig; P<0.01). The results of evolutionary convergence analysis also identified three instances of convergence on derived benthic-feeding head shapes (Figs 1 and 3).

The filmed species Amphiprion frenatus, Amphiprion ocellaris, Dascyllus aruanus, Chrysiptera cyanea, Pomacentrus moluccensis and Pomachromis richardsoni have PC1 scores that are intermediate between those of the two derived clusters in the 50 species head-shape space (Fig. 3A). The PC1 scores of planktivores and herbivores do not overlap for any of the 50 species examined, while those of omnivores, as might be expected, overlap with both of these groups. Although the derived planktivores constitute a group of species with distinct head shapes, no exclusive groups of either omnivores or herbivores can be distinguished (Fig. 3A).

As with the family-wide shape analyses, shape analyses of the heads, premaxillae and mandibles of the 10 species filmed during feeding distinguished between planktivores and fishes that that feed extensively from the benthos (herbivores and omnivores; Fig. 4). The head shapes of the species classified as planktivores in all evolutionary models (Table 1) were distinct from those of all other species (Fig. 4A). The head shapes of herbivores and omnivores (hereafter referred to collectively as ‘benthic feeders’) in evolutionary model 3.2 (Table 1), which was the best-supported evolutionary model (see below), were also distinct from each other (Fig. 4A), but this pattern did not hold true for premaxillary and mandibular shapes (Fig. 4D,G).

The premaxillary shapes of planktivores were distinct from benthic feeders (Fig. 4D). Except for S. nigricans, planktivores were distinguished from all other species by their PC1 scores alone (Fig. 4D). This axis accounted for a large percentage of the total premaxillary shape variation (82.75%) and was strongly associated with the relative lengths of the ascending and dentigerous arms (Fig. 4D,E). For all planktivores, the anterior and posterior extent of the ascending arms were similar in length, but in the case of S. nigricans the posterior side of the ascending arm of the premaxilla was much shorter than the anterior side (Fig. 4F). Long ascending arms permit greater upper jaw protrusion during feeding and this was evident in all planktivores (see below), but S. nigricans exhibited less upper jaw protrusion relative to its body length than did D. aruanus (Fig. S1), which has the shortest ascending arm length of all the species examined (Fig. 4F). Jaw manipulations of dissected, unfixed damselfishes indicated that the dorso-posterior side of the ascending arms of the premaxillae occlude with the maxillae when the mouth is completely closed. When the upper jaws were pulled forward, the extent of jaw protrusion was reached just before the tips of the ascending arms of the premaxilla became disarticulated from wing fossae of the articulation heads of the maxillae. This indicated that the length of the ascending arm measured along its posterior edge was directly related to maximum jaw protrusion distance and this ascending arm measurement was used in all statistical analyses. Planktivore mandible shapes were also distinct from those of benthic feeders (Fig. 4G), but as with head shape, this difference was not extensively due to shape differences associated with PC1 (Fig. 4A,G).

All individuals of every species used suction to capture food pellets, as movement of the pellets toward the fishes was observed coincident with mouth opening in all feeding events. In order to facilitate reading kinematic profiles, we grouped species by trophic niche (Fig. 5) in correspondence with the best-supported evolutionary model (model 3.2, see Table 1 and below). For kinematic profiles that depict all species, see Fig. S1.

Jaw movements of planktivores during feeding (upper jaw protrusion, gape and gape angle) strongly distinguished them from benthic feeders, and planktivores exhibited higher variation in these variables during the latter part of the feeding strike (Fig. 5A–C). Planktivore upper jaw protrusion was significantly greater immediately before and continuously after time zero (the time point at which a food pellet entered the mouth) in comparison to both omnivores and herbivores. Protrusion distance was much more variable among planktivores during this time period (Fig. 5A). Omnivores and herbivores had extremely similar protrusion profiles (Fig. 5A). Both planktivores and omnivores exhibited peak jaw protrusion at time zero, but the peak for omnivores was very small.

Planktivore gape profiles were largely different from those of the other species at time zero and afterwards (Fig. 5B), and planktivores showed the highest variability in gape during this time period. Gape differences during the latter part of the strike were only significant between planktivores and herbivores (Fig. 5B). All trophic groups exhibited peak gapes at time zero, and omnivore and herbivore gapes were significantly different at time zero. Omnivores exhibited extremely little variation in gape at this time point (Fig. 5B).

Planktivore gape angles were significantly greater than those of benthic feeders at time zero and afterwards. As with upper jaw protrusion and gape, planktivores showed the highest variability (Fig. 5C). Omnivores had significantly larger gape angles than herbivores at time zero. All trophic groups exhibited peak gape angles at time zero (Fig. 5C).

Body movements did not clearly distinguish between damselfishes in different trophic groups (Fig. 5D,E). Omnivores approached the food pellets at higher velocities than other fishes, but this difference was not significant in regard to planktivores for much of the strike (Fig. 5D). All species decreased body velocity immediately before food contact (Fig. 5D). Changes in cranial angle clearly did not differentiate among trophic groups. All groups showed an increase in cranial angle that began before time zero and peaked after time zero (Fig. 5E,F). Changes in cranial angle were small for all fishes (≤7 deg; Fig. 5F, Fig. S1). Our sampling rate was likely not high enough to accurately estimate acceleration and we omit characterizations or analyses of acceleration from all but our supplementary material (Figs S1–S3).

Contour map phylogenies that trace the evolution of the 18 variables listed in Table 2 (and acceleration as well) are depicted in Fig. S2. Table 3 shows which of these variables have undergone significantly correlated evolution. Evolutionary correlations plots are depicted in Fig. S3. The evolution of every variable was correlated with that of at least two other variables, but except for MV (eight significant correlations), those that exhibited the highest number of significant correlations (7–11) were all of those associated with upper and lower jaw shape, the movement of these structures during mouth opening and the integration of these movements (MJP, MG, MPR, APR, JA, LAA, LDA, ADA, OMA, FINT; Table 3). The large majority of the significant evolutionary correlations associated with these variables were with each other (Table 3). These findings strongly support our third hypothesis that jaw protrusion (MJP – maximum jaw protrusion) has evolved in correlation with a suite of additional morphological and functional traits associated with feeding.

The best-supported evolutionary model was the three-peak OU 3.2 model (see Table 1 for corresponding trophic classifications of species), which had an AICc score of −137.92 (Table 4). This finding supports our fourth hypothesis that the diversification of the form and function of damselfish feeding is best described by an evolutionary model of adaptation to three feeding niches: planktivory (pelagic feeding), omnivory (bentho-pelagic feeding) and herbivory (benthic feeding). The four-peak OU 4 model had the next strongest support (AICc=−135.30). This model differed from OU 3.2 in that D. aruanus was classified as a planktivore and S. nigricans was the sole representative of a fourth trophic group that was not used in any other model (detritivore). This AICc score difference of 2.62 indicates that although OU 3.2 is the preferred model, model OU 4 has support that is nearly as strong. None of the remaining models were well supported relative to these two. Three- and four-peak OU models had stronger support than the two-peak OU model, which in turn had better support than the single-peak OU model. The BM and EB models had the worst fit to the data. The evolutionary model-fitting analysis was robust to low sample size and exhibited high power to detect differences among models [null log-likelihood (LL)=−64.86, alternative LL=8.99; LL ratio=147.68; critical test statistic 101.40; power 0.97, P=0.008].

The trophic classifications used in OU 3.2 were also used in pANOVA testing. There were significant differences among trophic groups for eight of these variables (MJP, MG, MPR, APR, LAA, LDA, ADA, FINT; Table 5), and post hoc pANOVA results indicate that planktivores differed significantly from both omnivores and herbivores in their values for all of them except ADA (Table 5). Herbivores and omnivores exhibited no significant differences for any variable (Table 5). The finding that predominantly pelagic damselfishes (planktivores) have significant differences in jaw protrusion ability relative to species that feed extensively from the benthos (omnivores and herbivores) supports hypothesis 1. These results also support hypothesis 2 in that species capable of extensive jaw protrusion (planktivores) exhibited the highest levels of functional integration between upper and lower jaw movement (FINT). All of the variables that significantly distinguished between planktivores and benthic feeders were among those that have undergone correlated evolution with a large number of other variables (≥7) and most of their correlations were with each other (Tables 3, 5).

Analysis of our data supported all four of our hypotheses. We found that bentho-pelagic ecological divergence has been tightly linked to the evolution of jaw protrusion in damselfishes (hypothesis 1), that protrusion has evolved in correlation with multiple aspects of their trophic morphology (hypothesis 2), including functional integration between the upper and lower jaws (hypothesis 3), and that the evolution of their feeding mechanics is best described as adaptation to three trophic niches located at different points along the bentho-pelagic axis (hypothesis 4).

All aspects of head shape described by PC1 of the 50 species dataset (Fig. 3) are functionally associated with being able to feed from the water column versus the benthos. Terminal mouths allow for jaw protrusion in the direction of travel and aid with targeting free-swimming pelagic prey, while subterminal mouths allow feeding from the benthos while keeping the eyes positioned upward to detect predators (Konow and Bellwood, 2011). Larger jaw muscles aid the removal of tough algae from the substrate, while more robust jaw bones resist damage during hard biting. Insertion of the pars malaris near the maxillary–palatine joint facilitates rapid jaw closing while capturing elusive zooplankton, while an insertion farther from this joint provides a longer lever arm and facilitates hard biting. See Cooper and Westneat (2009) for further details.

Planktivores were significantly different from benthic feeders in regard to the extent and speed of upper jaw protrusion, the shape of the premaxilla, their ability to generate large gapes and the functional integration of mandible abduction with upper jaw protrusion during mouth opening (hypothesis 1; Table 5). These findings agree with modeling predictions of fish suction feeding mechanics which indicate that zooplankton capture is significantly affected by gape, mouth displacement speed and jaw protrusion distance (Holzman et al., 2012). Our data also support our second hypothesis that jaw protrusion ability (MJP) has evolved in correlation with a suite of additional morphological and functional traits associated with feeding, as the measures of shape and cranial movement that significantly differentiate benthic and pelagic feeders, including jaw protrusion, have largely evolved in correlation with each other (Tables 3, 5). Because the functional integration between lower jaw abduction and upper jaw protrusion (FINT) was one of the variables that has undergone significant, positively correlated evolution with MJP, these findings also support our third hypothesis that damselfishes capable of extensive protrusion will exhibit higher levels of functional integration between upper and lower jaw movement. Although Oufiero et al. (2012) found a significant correlation between jaw protrusion distance and attack speed in a kinematic study of serranid fishes, we did not find this to be true for the pomacentrids we examined. This may be due to the fact that the serranids were filmed while feeding on fishes, which is trophic specialization not seen among pomacentrids. Both of these studies were consistent in finding no correlation between the speed of jaw protrusion and attack speed (Table 3).

Although planktivores are distinct from benthic feeders in regard to their feeding biomechanics (Figs 3–5), the pANOVA results returned no significant differences between omnivores and herbivores for any of the variables considered individually (Table 5). However, when multiple variables were analysed together the evolutionary model that best fit the data was the three-peak OU 3.2 model (Table 4). This finding supports our fourth hypothesis that patterns in the diversification of damselfish feeding mechanics are best described by evolutionary models of adaptation to planktivory (pelagic feeding), omnivory (benthopelagic feeding) and herbivory (benthic feeding).

That planktivores differ markedly in both form and function from benthic feeders suggests that there may be trade-offs in the trophic biomechanics of damselfishes. In comparison to planktivores, benthic feeders not only have a reduced ability to protrude their jaws, but they also have larger biting muscles attached to more robust jaw bones in a manner that confers a higher mechanical advantage (Figs 2–4). Benthic-feeding damselfishes consume considerable amounts of tough algae and are better suited to produce strong bites than are planktivores. Planktivorous damselfishes must use high-speed, suction-producing strikes to capture elusive copepods whose long antennae are highly sensitive to the water displacement created by approaching predators (Day et al., 2015; Kiorboe and Visser, 1999; Webster et al., 2015; Yen et al., 2015). Their enhanced jaw protrusion distance and jaw protrusion speed facilitates this strategy (Holzman et al., 2008b; Motta, 1984).

Damselfishes in different trophic guilds were not distinguished by their location on the RSI (Table 5). As no damselfishes specialize on using body ram to capture elusive prey larger than zooplankton, the RSI may have little meaning for explaining niche-associated differences in their feeding mechanics. We specifically examined body ram in the present study, but Coughlin and Strickler (1990) described the rapid and extensive jaw protrusion of Chromis viridis as ram feeding. The ‘jaw ram’ that they observed is distinct from the movements measured here, and extensive jaw protrusion has been shown to enhance suction-dominant feeding (Ferry-Graham et al., 2001; Wainwright et al., 2001). Although it has been suggested that fishes will exhibit inverse abilities to utilize suction or ram to capture pelagic food items (Norton and Brainerd, 1993), it has been increasingly recognized that these abilities may not lie at either ends of a biomechanical spectrum (Ferry et al., 2015; Wainwright et al., 2001). Trade-offs for ram versus suction feeding may therefore not characterize feeding in many species.

Fishes that specialize on larger elusive prey (e.g. other fishes, shrimps) typically have large mouths and utilize fast ram strikes (Higham, 2007; Higham et al., 2007), but large mouths reduce suction feeding ability, which is only efficient within roughly one mouth diameter from the prey (Day et al., 2005; Wainwright et al., 2001). Although planktivorous damselfishes do have the largest standardized gapes (Table 5), their mouth sizes are nowhere near as large as those of piscivorous groupers (Serranidae), sunfishes (Centrarchidae) or cichlids (Cichlidae), and remain small enough to promote suction production (Cooper et al., 2010; Oufiero et al., 2012; Smith et al., 2015). We see no body ram specialists among the damselfishes and it has recently been suggested that the structure of their pharyngeal jaws may constrain their ability to specialize in larger prey, such as elusive fishes, whose capture is facilitated by body ram feeding strikes (McGee et al., 2015a).

Exclusively morphological studies have shown that damselfish trophic evolution can be characterized as repeated adaptation to planktivory, omnivory and herbivory (Cooper and Westneat, 2009; Frédérich et al., 2016, 2013; Olivier et al., 2014) and this description is further supported by our combined analyses of both morphological and kinematic data (Table 4). Heavy feeding on detritus has evolved at least three times in herbivorous damselfishes (Stegastes nigricans, Hemiglyphidodon plagiometopon and Dischistodus perspicilliatus; Wilson and Bellwood, 1997) and the four-peak OU model, which included detritivory, had a fit to our data that was nearly as good as OU 3.2 (Table 4). This may, however, simply represent the fact that the trophic morphology of S. nigricans, our only detritivorous species, is distinct from those of the other herbivores examined here (Fig. 3). All diet studies of this species have reported extensive feeding on algae (Table 1) (Letourneur et al., 1997; Sano et al., 1984; Wilson and Bellwood, 1997).

Premaxillary shape is of strong importance to the evolution of jaw protrusion (Table 5) (Hulsey et al., 2010) and this structure represents a potential target for evo-devo studies of the developmental determinants of protrusion ability. Developmental studies have expanded to include many cypriniform and acanthomorph fish species in recent years (e.g. zebrafish, medaka, barbs, halfbeaks, flounder, cichlids, damselfishes, etc.; Cooper et al., 2013; Gunter et al., 2014; Kavanagh and Alford, 2003; Langille and Hall, 1987; Le Pabic et al., 2016; Power et al., 2001; Shkil and Smirnov, 2015). Therefore, it may be possible to investigate whether the convergent evolution of jaw protrusion has involved similar developmental changes in these two lineages.

Our observations that: (1) species with highly protrusile jaws have significantly higher levels of functional integration between the upper and lower jaws (Table 5) and (2) functional integration of the jaws has undergone correlated evolution with all aspects of jaw form and function that distinguish between planktivores and benthic feeders (Table 3) demonstrate the importance of the integration of upper and lower jaw movement to the trophic divergence of the Pomacentridae. Whether fishes with highly protrusile jaws have higher levels of covariation between upper and lower jaw shape (i.e. greater morphological integration) is a testable hypothesis. Support for this hypothesis would indicate that increased protrusion would require the evolution of developmental mechanisms that impose higher levels of upper and lower jaw shape covariation (i.e. developmental integration; Hallgrímsson et al., 2009).

Maximum jaw protrusion distances strongly distinguish between planktivorous and benthic-feeding damselfishes. Multiple components of pomacentrid bite mechanics have evolved in correlation with protrusion ability, including functional integration between the upper and lower jaws. The evolution of their trophic morphology corresponds to a pattern of adaptation to three feeding niches located at different points along the bentho-pelagic axis. The significant differences in both form and function that distinguish planktivorous damselfishes from benthic feeders suggest that trade-offs in regard to bite speed and bite force may have played an important role in the evolution of pomacentrid trophic mechanics. Finally, we suggest that the developmental controls of premaxillary shape constitute an important target for evo-devo studies of fish feeding.

We gratefully acknowledge the Lizard Island Research Station for their support. The Australian Museum assisted with shipping specimens to the United States. Holly Milewski assisted with analysing video data for Amphiprion ocellaris.