Phylogenetic evidence for recent diversification of obligate coral-dwelling gobies compared with their host corals.

The rich diversity of coral reef organisms is supported, at least in part, by the diversity of coral reef habitat. Some of the most habitat specialised fishes on coral reefs are obligate coral-dwelling gobies of the genus Gobiodon that inhabit a range of coral species, mostly of the genus Acropora. However, the role of this specialised pattern of habitat use in the evolution of coral-dwelling gobies is not well understood. Diversification of coral-dwelling gobies may be driven by the diversification of their host corals (cospeciation), or alternatively, diversification of these fishes may have occurred independently of the diversification of host corals. The cospeciation hypothesis assumes similar timing in evolution of the gobies and their host corals. We used four genes for each group and the available fossil records to reconstruct and date phylogenies for 20 species of Gobiodon from the Indo-Pacific and the Red Sea, and for 28 species of the coral genus Acropora. Our results indicate that Gobiodon diversified mostly in the last ∼5My, whereas Acropora corals have consistently diversified since the Eocene, making the hypothesis of cospeciation untenable. The fully resolved molecular phylogeny of the genus Gobiodon is in part at odds with previous analyses incorporating morphological data and indicates that some morphological traits form paraphyletic clades within Gobiodon. Our phylogeny supports a hypothesis in which Gobiodon diversified in the Indo-Pacific Ocean and then radiated recently, with multiple new variants found in the Red Sea.

Introduction

Species interactions can influence evolution and result in coevolved systems (Thompson, 2009). If interactions between species are close enough, the organisms involved may speciate at similar evolutionary times, so a reconstruction of their evolutionary histories would show congruent events of speciation; a pattern known as cospeciation (Paterson and Banks, 2001). Phylogenetic inferences have been used to study cospeciation in closely interacting groups of organisms, mainly from parasitic and mutualistic associations (e.g. fish and their parasites, Huyse and Volckaert, 2005; figs and fig wasps, Weiblen and Bush, 2002; salamanders and their viruses, Storfer et al., 2007; echinoderms and annelids, Lanterbecq et al., 2010; yucca plants and yucca moths, Althoff et al., 2011). Although evidence for cospeciation events is most often supported, the assumption of similar timing of evolution of interacting groups is rarely tested (Paterson and Banks, 2001). Current techniques to date evolutionary events provide a tool to test this assumption robustly and may give better insight into the prevalence of cospeciation.

Interactions between reef fishes and corals are important to reef ecosystems both ecologically and evolutionarily (Jones et al., 2004; Graham et al., 2006; Cowman and Bellwood, 2011), therefore it is important to understand the history of their associations (Bellwood and Wainwright, 2002; Rocha and Bowen, 2008; Kiessling et al., 2010). Reef fishes and reef building corals are both known for their extraordinary diversity; however, the evolutionary links between these two groups of coral reef organisms is not fully understood. Some of the most habitat specialized fishes on coral reefs are from the genus Gobiodon, which are obligate coral-dwelling gobies that mostly inhabit coral colonies from the genus Acropora (Munday et al., 1997, 1999). Gobiodon species are highly selective among coral species (Munday et al., 1997; Munday, 2004a; Dirnwoeber and Herler, 2007) and they compete strongly for access to preferred coral species, both within and between species (Munday et al., 2001; Hobbs and Munday, 2004). The fitness of coral gobies depends on the availability of coral host species (Munday, 2001; Caley and Munday, 2003; Herler et al., 2011), which indicates that the association has evolutionary ramifications for the fish. Given that the genus comprises at least 20 species, many of which remain faithful to their preferred coral host across reef systems (Munday, 2000, 2002), their evolution and interaction with Acropora provides an intriguing case study of reef fish evolution.

Ecological aspects of the association between Gobiodon species and Acropora corals may suggest a linked evolutionary history. Coral-dwelling gobies depend on coral colonies as a source of shelter, food and breeding sites (Munday et al., 2001; Hobbs and Munday, 2004; Brooker et al., 2010). Most species of Gobiodon have specific preferences and inhabit a limited number of the Acropora coral species present on coral reefs (Munday et al., 1999, 2004). Furthermore, habitat use of some Gobiodon species varies little across geographic regions in the Indo-West-Pacific (Munday, 2002), indicating that habitat choice can be a highly constrained trait. Gobiodon also provide benefits to their host corals; for instance, chemical cues released by Acropora corals that have come into contact with toxic algae trigger the gobies to feed on the algae and prevent it overgrowing the coral (Dickson and Hay, 2012). Similarly, gobies may protect the coral from predation by butterflyfishes and other corallivorous fishes (Dirnwoeber and Herler, 2012). This close association between Gobiodon and Acropora corals is suggestive of a linked evolutionary history.

Past phylogenetic analyses of the genus Gobiodon have been based on morphological and mitochondrial genetic data (12S and 16S rRNA; Harold et al., 2008). The genus is monophyletic (Harold et al., 2008; Herler et al., 2009), although previous phylogenetic analyses using maximum parsimony resulted in low resolution for the internal nodes (Harold et al., 2008). Morphological features like body shape and osteological structures have been used in attempts to improve this resolution (Harold et al., 2008). Nevertheless, there remains considerable uncertainty about the phylogenetic relationships among Gobiodon species and there is some discordance between morphological and molecular results (Harold et al., 2008). The weak resolution of phylogenetic relationships in previous analyses means that they have limited power for testing evolutionary hypotheses about diversification in this group of fishes. The timing of diversification of the genus Gobiodon has been estimated in a previous study, suggesting that the group started diversifying around 10 Mya (Herler et al., 2009). However, this study used a subset of just eight species, used only one mitochondrial gene, and assumed a molecular clock, precluding a robust comparison with the timing of diversification in their host corals from the genus Acropora.

The diversification of scleratinian corals extends to the last 200 My (Simpson et al., 2011). To date, the timing of evolutionary events in the genus Acropora has been hypothesized utilizing phylogenies based on cladistic and maximum-likelihood analyses (Wallace, 1999; van Oppen et al., 2001). Most of the diversification of the Acropora genus was thought to have occurred in the Pliocene and Pleistocene (Wallace, 1999; Van Oppen et al., 2001); however, recent fossil findings suggest that the Cervicornis species group was already present in the Lutetian (∼45 Mya; Wallace and Rosen, 2006; Wallace, 2008). Fossils from the Hyacinthus and the Aspera species groups also suggest an earlier divergence of the genus (Wallace, 2008). The rate of molecular evolution of Acropora corals is slow and there are mechanisms in place that retard it (Van Oppen et al., 1999). These mechanisms have been suggested to cause unusually extended longevity, and include asexual reproduction, slow growth and the lack of a mortal soma (Hellberg, 2006). Similarly, Acropora corals are likely to undergo hybridization and introgression, which may also lead to a slower rate of molecular substitutions and reduce the rate of extinction (Willis et al., 2006; Richards et al., 2008).

Studies on the evolutionary history of reef fishes have suggested that the interaction between fishes and coral reefs became common soon after the Cretaceous-Tertiary (K/Pg) boundary (Cowman and Bellwood, 2011). However, little is known about the evolutionary history of the interaction between gobies and corals. The present study undertakes an evolutionary analysis of the association between Gobiodon and Acropora by independently estimating the dates of diversification for both groups. Specifically, we test if the diversification times of Gobiodon and Acropora overlap, which is a necessary assumption of the cospeciation hypothesis that is rarely tested (Paterson and Banks, 2001). We constructed the most complete genetic dataset for Gobiodon to date, using mitochondrial and nuclear markers for 20 species from the Indo-Pacific Ocean and the Red Sea. Similarly, we compiled data for four markers of Acropora of mitochondrial and nuclear origin, including 12 of the species most commonly inhabited by Gobiodon. Using available fossil data we inferred a phylogeny of Gobiodon and estimate the timing of first appearance of Gobiodon and Acropora to test for co-incidental dates of evolution. To gain further insight into the evolutionary dynamics of the two groups, we use a Bayesian approach to test the plausibility of evolutionary models of pure speciation (Yule process) against a model including both speciation and extinction (birth–death process). Comparing the dates and dynamics of evolution of these two groups is instrumental for future studies of the diversification of Gobiodon and the origins of their association with Acropora.

Materials and methods

Taxon sampling

Our phylogenetic analysis comprised twenty Gobiodon species (Harold et al., 2008), including 6 recognized, but as yet undescribed species (G. sp. A, B, C, D; Munday et al., 2004, 1999; and G. sp. 1 and 2; Herler et al., 2009). Additionally, we included samples from the Red Sea of three species that are widespread in the Pacific Ocean (G. rivulatus, G. histrio and G. citrinus; Harold et al., 2008; Munday et al., 1999; Herler et al., 2009). Amblyeleotris sp., Ctenogobiops sp., and Paragobiodon xanthosomus were used as outgroup species (Table 1). Specimens of Indo-Pacific species were collected by PLM at Lizard Island on the Great Barrier Reef, Australia, and in Kimbe Bay, Papua New Guinea (Table 1). Specimens of Red Sea species were collected by JH from the Gulf of Aqaba (Dahab, Egypt), the northern Red Sea (Marsa Alam, Egypt), and the southern Red Sea (Dahlak Archipelago, Eritrea). Specimens were collected following anaesthetization with dilute solution of clove oil (Munday and Wilson, 1997) and preserved in 80% alcohol prior to molecular analyses.

Table 1

Gobiodon and outgroup species collected and examined in this study indicating sample locations and references. PO-Pacific Ocean, RS-Red Sea, GBR-Great Barrier Reef, PNG-Papua New Guinea.

Taxon	Location	Reference
Gobiodon acicularis	Lizard Island, GBR, Australia	Munday et al. (1999)
G. axillaris	Lizard Island, GBR, Australia	Munday et al. (1999)
G. brochus	Lizard Island, GBR, Australia	Munday et al. (1999, 2004)
G. ceramensis	Lizard Island, GBR, Australia	Munday et al. (1999)
G. citrinus PO	One Tree Island, GBR, Australia	Munday et al. (1999)
G. citrinus RS	Northern Red Sea (Gulf of Aqaba, Marsa Alam), southern Red Sea (Dahlak Archipelago)	Herler et al. (2009)
G. erythrospilus	Lizard Island, GBR, Australia	Munday et al. (1999, 2004)
G. histrio PO	Lizard Island, GBR, Australia	Munday et al. (1999, 2004)
G. histrio RS	Northern Red Sea (Gulf of Aqaba, Marsa Alam)	Herler et al. (2009)
G. oculolineatus	One Tree Island, GBR, Australia	Munday et al. (1999)
G. okinawae	One Tree Island, GBR, Australia	Munday et al. (1999)
G. quinquestrigatus	Lizard Island, GBR, Australia	Munday et al. (1999)
G. reticulatus RS	Northern Red Sea (Gulf of Aqaba, Marsa Alam), southern Red sea (Massawa)	Herler et al. (2009)
G. rivulatus PO	One Tree Island, GBR, Australia	Munday et al. (1999)
G. rivulatus RS	Northern Red Sea (Gulf of Aqaba, Marsa Alam), southern Red sea (Massawa)	Herler et al. (2009)
G. spilophthalmus	Lizard Island, GBR, Australia	Munday et al. (1999)
G. unicolor	One Tree Island, GBR, Australia	Munday et al. (1999)
G. sp. A	Lizard Island, GBR, Australia	Munday et al. (1999)
G. sp. B	Bootless Bay, PNG	Munday et al. (1999)
G. sp. C	Kimbe Bay, New Britain, PNG	Munday et al. (1999)
G. sp. D	Kimbe Bay, New Britain, PNG	Munday et al. (1999)
G. sp. 1 RS	Northern Red Sea (Gulf of Aqaba), southern Red Sea (Dahlak)	Herler et al. (2009)
G. sp. 2 RS	Northern Red Sea (Gulf of Aqaba, Marsa Alam)	Herler et al. (2009)
Outgroup taxon
Paragobiodon xanthosomus	Lizard Island, GBR, Australia
Amblyeleotris sp.	Lizard Island, GBR, Australia
Ctenogobius sp.	Lizard Island, GBR, Australia

Laboratory procedures

Total DNA was extracted from Gobiodon tissues using standard salt-chloroform and proteinase K digestion extraction procedures (Sambrook and Russell, 2001). Four loci, three mitochondrial (12SrRNA, 16SrRNA, cytochrome b) and one nuclear, S7 ribosomal protein gene Intron1 (S7I1), a gene required for assembling RNA (Chow and Hazama, 1998; Maguire and Zimmermann, 2001), were sequenced. In addition to the four markers utilized in this study, we designed specific Gobiodon primer sequences for cytochrome b (Supporting Information Table 1). The Gobiodon primers were used interchangeably with published cyt b primers either replacing L14841 or H15149 (Kocher et al., 1989) in PCR reactions. Each 20 μl PCR reaction volume contained 2.5 mM Tris–Cl (pH 8.7), 5 mM KCl(NH4)2SO4, 200 μM each dNTP, MgCl2 ranging from 1.5 mM to 4 mM, 10 μM each primer, 1 unit of Taq Polymerase (Qiagen) and 10 ng template DNA. Amplifications followed the same basic cycling protocol: an initial denaturing step of 2 min at 94 °C, followed by 35 cycles, with the first 5 cycles at 94 °C for 30 s, 30 s at primer specific annealing temperatures (T; see SI Table 1); followed by 1 min 30 s extensions at 72 °C and the remaining 30 cycles were performed as before, but at T −2 °C. PCR products were purified by isopropanol precipitation and sent to Macrogen Inc. (Korea) for sequencing on an ABI 310 XL sequencer using ABI dye-terminator chemistry.

Sequence data compilation

Three separate data sets were compiled for analyses: one for the Gobiodon phylogeny, a second and larger dataset, hereafter called the Gobiiformes with representative species of the families Gobiidae, Gobionellidae, and Apogonidae, which was required to estimate the age of the emergence of Gobiodon and a third dataset for estimation of the diversification age of Acropora corals. The latter two datasets were downloaded entirely from GenBank (see SI Table 2).

Gobiodon and Gobiiformes

Two specimens were sequenced for most species of Gobiodon. An additional tip in the phylogeny and pair of specimens was included for gobies with representatives in the Indo-Pacific and the Red Sea. A single specimen was used for G. ceramensis, G. citrinus (from the Indo-Pacific), and G. spilophthalmus. The consensus of duplicate sequences was created using Sequencher 4.5 (Gene codes corporation). The resulting dataset was automatically aligned using ClustalX (Thompson et al., 1997) and corrected manually using Se-Al version 2.0 available at http://evolve.zoo.ox.ac.uk (Rambaut, 1996). Sequences have been deposited at GenBank accession numbers KC894468-KC894517. Additional sequences for most of 12S- and 16SrRNA, and cyt b of six Red Sea species or variant populations of Gobiodon were downloaded from GenBank Accession No.’s: 12SrRNA: EF540558-EF540584, FJ617027-FJ617038, and FJ617041-FJ617046;16SrRNA: EF443263-EF443264, EF443267-EF443268, EF463067-EF463076, EF527238-EF527252, EF527254, FJ617067-FJ617078, and FJ617081-FJ617086; cyt b: FJ617107-FJ617118, and FJ617121-FJ617126 (Harold et al., 2008; Herler et al., 2009). Three genes, the mitochondrial 12S- and 16SrRNA and the nuclear intron (S7I1) were partitioned into putative stem (conserved 12S = 280 bp, 16S = 396 bp, S7I1 = 540 bp) and loop (hypervariable 12S = 104 bp, 16S = 126 bp, S7I1 = 131 bp) regions; the fourth gene cyt b was partitioned into 1st and 2nd codon positions combined as conserved (247 bp) and 3rd codon positions as variable regions (123 bp). This biologically realistic partition scheme can reasonably capture the heterogeneity in the data (Brandley et al., 2005; Brown and Lemmon, 2007), and has proven effective in other studies of reef fishes (Bellwood et al., 2010; Choat et al., 2012). In total, eight separate gene partitions (g1–g8) were identified and each region was examined for its best fitting model using MrModeltest version 2.2 and Aikaike information criterion (AIC) (Nylander, 2004; Nylander et al., 2004). The eight separate gene partitions, each with their specific model, were subsequently concatenated for further phylogenetic analyses.

The larger Gobiiformes dataset contained sequences of four mitochondrial markers (CO1, ND1, ND2, and cyt b; Thacker, 2009) and was examined primarily for molecular age estimations (SI Table 2). Two of the resulting age estimates were used as secondary calibrations to estimate the age of the most recent common ancestor (MRCA) to the Gobiodon genus.

Acropora corals

Within the coral genus Acropora, we only used genes that were available for a wide range of species, which resulted in a dataset of two mitochondrial (control region and NAD5) and two nuclear markers (PaxC 46/47 intron and the Calmodulin CaM-encoding gene). The final dataset included 28 species from the genus Acropora, with four species from the genus Isopora as outgroup. Four individuals for each species were used for analyses of each gene to observe potential introgression patterns as seen in previous phylogenetic studies of this genus (e.g. Van Oppen et al., 2001 or Richards et al., 2010). We included representative species from 14 of the 19 recognized species groups (Wallace, 1999), which accounts for close to 17% of the current valid Acropora species (Veron, 2002) and represents about 60% of the coral species inhabited by Gobiodon (SI Table 2).

Phylogenetic inferences

Six gene-specific models of substitutions were chosen based on AIC and applied in the following phylogenetic analysis (see SI Table 3). The model selection for pMM Bayesian analysis only requires a general ‘form’ of the model, as the Markov chain integrates uncertainties of the parameter values (Nylander et al., 2004). Therefore, the base frequency was set to fixed = equal for three of the eight gene partitions (12 and, 16S stems, S7 loop), while the remaining five partitions (12 and 16S loops, cyt b conserved and variable and S7 stem) were set with unequal base frequencies.

Bayesian inference (BI) phylogenetic analyses were implemented in Mr. Bayes version 3.1.2 (Huelsenbeck and Ronquist, 2001) using CIPRES Portals (Miller et al., 2009), accessed at the following URL site http://www.phylo.org/sub_sections/portal. The analysis of the concatenated data used a partition-mixed model method (pMM) utilising the identified locus-specific substitution models. Five Bayesian pMM analyses were performed using Markov chain Monte Carlo (MCMC) simulations with four chains of 5,000,000 generations each, sampling trees every 500 generations. Appropriate mixing was reached after 20,000 generations, visualised in Tracer v1.5 (Rambaut and Drummond, 2007), and a 50% majority-rule consensus tree was computed using the best 500 post burn-in trees from each run. Three outgroup species, Amblyeleotris sp. and Ctenogobiops sp., and Paragobiodon xanthosomus were used to root resulting trees. The single best tree was selected for molecular dating.

Maximum parsimony (MP) analyses were implemented in PAUP* 4.0b10 (Swofford, 1998) using heuristic search methods with 1000 pseudo-replicate bootstraps, tree-bisection-reconnection branch swapping and random addition of taxa. Two separate heuristic MP runs were performed. In the first all sites were treated equally, and in the second sites were weighted 2:1 according to gene partitions; sequences from the mtDNA (12S, 16S and cyt b; 1276 bp) were given a weight of 1, and the nuclear gene S7I1 (671 bp) was given a weight of 2 in order to compensate for the smaller number of characters in the single nuclear gene used. A 50% majority-rule consensus tree was generated from all shortest trees obtained.

Sequences acquired from GenBank (SI Table 2) for datasets of Gobiiformes and Acropora were aligned using the Muscle algorithm (Edgar, 2004) and then manually checked with the software SeaView v4.3 (Gouy et al., 2010). This was followed by a process of substitution model selection with both the Akaike Information Criterion and the Bayesian Information Criterion in the software JModelTest (Posada, 2008; SI Table 3).

The initial phylogeny for the Gobiiformes dataset was estimated using BEAST v1.7.2 with 10M steps in a single chain, which recovered the same species relationships as the original study by Thacker (2009). The topology with the higher product of posterior probabilities from this run was used as the initial tree to estimate ages of Gobiiformes.

Acropora

We assessed phylogenetic relationships of each gene separately to determine the best partitioning model. Five maximum likelihood runs were performed for each gene with one hundred replicates using Garli v1.0 (Zwickl, 2006). The best tree for each gene was chosen using the Shimodaira–Hasegawa test as implemented in the Phangorn v1.4-1 package in R (Shimodaira and Hasegawa, 2001; Schliep, 2011; www.r-project.org; R Development Core Team 2012), and the congruence of trees between genes was used to determine the best partitioning scheme. Using these partitions, we estimated separate phylogenies with a 20M step Bayesian analysis in the program MrBayes v3.0 (Huelsenbeck and Ronquist, 2001). As with the Gobiiformes the final topologies were used as initial values to estimate dates in the phylogeny of Acropora corals.

Molecular dating

Gobiiformes and Gobiodon

Relaxed molecular clock models were applied to all datasets for diversification estimates while accounting for associated uncertainties (Graur and Martin, 2004). Fossil and secondary calibrations were given a lognormal distribution, which is more appropriate than other forms of priors for these two types of calibrations (Ho and Phillips, 2009). Given the fossil record for Gobiidae is limited, we used secondary calibrations acquired from our dating estimates of the Gobiiformes phylogeny (SI Fig. 1). A fossil prior of the genus Pomatoschistus sp. (Carnevale et al., 2006) was used as a prior for the root of the suborder Gobioidei. Similarly, the fossil of Eosphaeramia sp. for apogonids (Bannikov, 2008) and the K/Pg boundary for the root were also used as date priors (SI Table 4). The K/Pg boundary event was placed at the tree root with a normal distribution to allow age estimates to capture both older and younger ages around the 65.5 My event. Finally, the phylogenetic and substitution model estimations from Section 2.5.2 were used as initial priors and models to estimate divergence dates using a Markov Chain Monte Carlo (MCMC) Bayesian analysis implemented in BEAST v1.7.2 (Drummond and Rambaut, 2007) with 20M steps.

The age estimates of interest to employ as secondary calibrations on the phylogeny of Gobiodon are 22.01 Mya (95% HPD = 15.65–28.77) for the divergence between Ctenogobiops sp. and Ambyeleotris sp., and 30.58 Mya (95% HPD = 25.14–35.88) at the root of these two outgroup genera with Gobiodon (SI Fig. 4). These two age estimates were set as priors for Gobiodon and evaluated in BEAST v1.7.2 with three 20M step runs including parameters for the analyses as outlined above.

The oldest coral fossil for the Acropora group III (sensu Van Oppen et al., 2001) was identified as Acropora wilsonae (∼44.4 Mya; Wallace, 2008). This calibration prior also comprised the oldest known fossil from the Cervicornis group, Acropora alvarezi, found at the same site as Acropora wilsonae (Wallace, 2008). In addition, the oldest known Acropora fossil has only been identified to the level of genus (Carbone et al., 1993), and was included as a calibration at the root of the genus (SI Table 4). Two runs were performed with 50M step chains each to reach satisfactory mixing of the distribution. Files of replicate runs were combined using LogCombiner v1.7.0 and maximum clade credibility trees were identified using TreeAnnotator v1.7.0 (Drummond and Rambaut, 2007).

Comparison of date estimates

The 95% highest posterior density intervals (HPDIs) for the divergence estimates were compared to evaluate the congruence in dates of cladogenesis between Gobiodon and Acropora. To identify the most appropriate tree construction prior, between a prior that includes extinction and one that excludes it, we used the bayes factors of the marginal likelihoods as estimated in the software Tracer v1.5 (Drummond and Rambaut, 2007). These two models (also called Yule process and birth–death respectively) provide an indication of the evolutionary dynamics taking place in each taxon. The proper mixing of bayesian runs were observed in Tracer v1.5, and only runs with effective sample sizes (ESS) above 200 for all parameters were examined.

Results

We examined 1947 bp of Gobiodon species sequences. Phylogenetic inference for this group showed congruence between topologies obtained from different analyses (Bayesian inference and maximum parsimony; Fig. 1). Gobiodon was monophyletic with total support and contained four distinct clades, three of which had Bayesian Inference and Maximum Parsimony support of ⩾87% (Fig. 1). Short branches at the base of all four clades suggest an accelerated radiation in this genus. Clade I contained two sister species groups containing five described species and the Red Sea variant of G. citrinus, which is more divergent from its Pacific Ocean counterpart than are the other sister species from each other in this clade. Clade II retrieved a strongly supported group of two sister species, G. sp. A and G. sp. B, and their sister, G. brochus, whilst the placement of G. sp. C was not bilaterally supported (Fig. 1). Clade III was strongly supported and had two sister groups, each including a species from the Red Sea or a Red Sea population of an Indo-Pacific species (G. sp. 2 with G. axillaris and G. unicolor in one group, and G. histrio from the Pacific and Red Sea in the other). Clade IV had very strong support throughout and contained two lineages with three sister species pairs. The first lineage contained G. oculolineatus and G. sp.1 from the Red Sea as its closest relative, which share a common ancestor with the Red Sea species G. reticulatus. The second lineage consisted of G. quinquestrigatus and G. sp. D as sister species, and G. rivulatus (Fig. 1).

Fig. 1

Phylogenetic inferences of the genus Gobiodon based on comprehensive taxon sampling with three known outgroup species, Paragobiodon xanthosomus, Amblyeleotris sp. and Ctenogobiops sp., obtained by Bayesian and maximum parsimony (MP) analyses of four loci (12S and 16SrRNA, cytochrome b, and nDNA S7I1). Topology of best Bayesian tree (consensus of 500 post burn-in trees from each run) with posterior probabilities (%) and bootstrap support (>50%) of MP (1000 bootstrap replicates) are indicated. Asterisks (*) represent 100% posterior probability/bootstrap support respectively. Four main clades are indicated. RS: Red Sea variant, PO: Pacific Ocean (western region) variant.

Chronological comparisons

Our re-analyses of Gobiiformes with 4362 bp, to derive age estimates, required the same substitution model (GTR + I + G) and retrieved the same phylogenetic topology as was retrieved by Thacker (2009; SI Figure 1). Our results show that the diversification of the present species of Gobiodon is recent, occurring mostly during the Pleistocene (Fig. 2). The youngest species are the sister species G. sp. A and G. sp. B, which diverged during the late Pleistocene, <1 My (95% HPDI ∼0.5–2 My). At the other end of the age spectrum, the ancestral species that gave rise to the present genus of Gobiodon emerged during the late Miocene, approximately 8 My (95% HPDI of ∼4.5–12 My) (Fig. 2). Only the oldest three Gobiodon species, depicted on the chronogram (G. sp. C, G. brochus and G. citrinus; Fig. 2), are of possible late Miocene origin (<10 Mya).

Fig. 2

Chronogram based on BEAST Markov chain Monte Carlo runs with 95% highest posterior density interval (HPDI) in million years (My). Numbers on branches represent posterior probabilities of the topology inferred separately by BEAST.

Nuclear and mitochondrial datasets of Acropora were analyzed separately to infer chronograms using Bayesian analyses (SI Fig. 2), due to their topological incongruence in exploratory Maximum Likelihood analyses. Despite this, 95% highest posterior density intervals (HPDIs) of date estimates of both datasets overlapped for all species except Acropora divaricata (Fig. 3) and ranged from a minimum of 5.8 to a maximum of 52.9 My. Eight of the twelve species analyzed that have association with Gobiodon had a minimum 95% HPDI of >10 My.

Fig. 3

MCMC Bayesian results for inferred dates of emergence of species of Gobiodon and coral hosts of Acropora. The bars indicate the 95% Highest Posterior Density Interval (HPDI) with the mean estimate indicated by an inner bar.

The origin of the genus Acropora was likely during the early Eocene and more recent than the K/Pg boundary (Fig. 3). The diversification of the genus for both mitochondrial and nuclear data appears to start soon after the K/Pg boundary with the separation of Clade I between 44 and 56 Mya (SI Fig. 2). Although the age estimates for the first appearance of species were generally more recent as presented by the mitochondrial data, there is consistent overlap between the two datasets. Another shared feature among Acropora datasets is the wide HPDIs, which is associated with the high phylogenetic uncertainty. Nevertheless, the ages for first emergence of Acropora species that have Gobiodon associated with them in the present was mostly during the late Eocene. The oldest mean age of diversification of Acropora that host Gobiodon in the present are A. gemmifera, A. cerealis and A. digitifera. While most of the mean ages for the emergence of Acropora species occur during the Oligocene and Eocene epochs, they are all likely to have emerged earlier than 5 Mya (Fig. 3). These results show that diversification between gobies and their coral associates has occurred at a fundamentally different time scale.

Diversification comparison

The posterior distribution for the yule. BirthRate parameter, which is proportional to the net speciation rate, was around twice the width for Gobiodon (0.08–0.22) compared with Acropora (0.06–0.097), indicating larger variability in the speciation rate in gobiid fishes than in the acroporid corals (Fig. 4). The comparison between tree construction priors that both included and excluded extinction processes (birth–death and Yule respectively) produced a bayes factors value of 31.58, in favor of the pure birth process in Acropora, suggesting that extinction has not been important in the evolutionary history of Acropora. In contrast, bayes factors produced from Gobiodon data generated a value of 0.68, suggesting that the birth–death model could not be rejected and that extinction has likely been important in the evolutionary history of Gobiodon. The likely presence of both speciation and extinction events in the evolutionary history of Gobiodon is consistent with the greater range of speciation rates evident in Gobiodon compared to that of Acropora.

Fig. 4

Inferred posterior distribution from Bayesian analysis of the parameter yule. Birthrate (proportional to the rate of speciation) for the genera Gobiodon and Acropora.

Discussion

We resolved the phylogenetic relationships within the genus Gobiodon and estimated ages of diversification for these species and the Acropora corals they inhabit. Most of the diversification of the genus Acropora occurred during the Oligocene (∼34–23 Mya), although the whole diversification of the species studied, including those associated with species of Gobiodon, extended from the mid-Eocene to the Miocene (∼49 to ∼5 Mya). In contrast, fishes of the genus Gobiodon, only started diversifying in the late Miocene (∼9 Mya) and most of their diversity is of Pleistocene origin (∼2.8–0.1 Mya). Therefore, the diversification of Gobiodon and Acropora mostly occurred at different times, indicating a lack of evolutionary concordance. Although these groups have a close ecological interaction, their diversity is likely to have independent evolutionary origins and the hypothesis of cospeciation between them is not supported.

The evolution of

Our analysis shows that the genus Gobiodon has two well-supported major clades, each with two minor clades. Our inferred Gobiodon phylogeny retrieved similar relationships between species at the tips of the tree as the molecular analysis of Harold et al. (2008), albeit with additional species to the earlier analysis. However, using four molecular markers (including one nuclear marker) we were able to resolve the relationship of the four clades to each other, which was not possible in the earlier analysis of Harold et al. (2008) that used just two mitochondrial markers.

Our molecular analysis suggests that sister species G. sp. A and G. sp. B are most closely related to G. brochus, all of which form a monophyletic clade with G. sp C. This contrasts with the molecular analysis of Harold et al. (2008) that was not able to resolve robustly the relationship between G. brochus, G. sp. C and the sister species G. sp A and G. sp. B and which, with the inclusion of morphological data, suggested that this group of species was not monophyletic. Our analysis provides strong support for a monophyletic clade.

Previous combined molecular and morphological analyses by Harold et al. (2008) indicated that G. erythrospilus and G. histrio were closely related to what we refer to as Clade II, which is not the case in our molecular analyses. This means that some of the morphological features incorporated in previous analyses such as pronounced deep body shapes (G. brochus, G. histrio, G. unicolor, G. erythrospilus, G. axillaris, G. sp. A, G. sp. B, and G. sp. C) and an interopercular-isthmus groove (same species as deep body shape, but excluding G. axillaris), form paraphyletic clades. The distribution of these morphological features among clades is consistent with a rapid evolution of Gobiodon and high morphological plasticity, and is also apparent in the unique features that separate this genus from other gobies, for example deep body shapes and a lack of scales (Herler et al., 2011). Deep body shape has evolved in some other coral-dwelling fishes that are not closely related to Gobiodon, such as crouchers (Caracanthidae) and damselfishes (Pomacentridae), indicating that this morphological trait is adaptive for fishes that live permanently among the branches of coral colonies and is not likely to be a reliable indicator of phylogenetic relatedness.

In our analysis, species and populations sampled from the Red Sea were sister to species and more widespread populations of Indo-Pacific Gobiodon, indicating that the Red Sea was repeatedly and independently colonized by Indo-Pacific Gobiodon species. This is not surprising considering that the Red Sea is one of 10 coral reef fish hotspots of endemism in the Indo-Pacific (Allen, 2008). It is also indicative that the genus may have diversified in the periphery of its range, as has been shown with other reef fish clades (e.g. Winters et al., 2010; Hodge et al., 2011). Whole clades of Gobiodon, like that of G. reticulatus, G. sp.1, and G. oculolineatus, may have originated in the Red Sea as the region was recently found to be a significant source of marine fish fauna (DiBattista et al., 2013). However, the genetic extent and age of colonization events can only be described with further phylogeographic and population-level analyses.

Neither of the tree construction models, pure birth (Yule) or birth–death process, was better suited for explaining the evolution of Gobiodon, suggesting that either extinction events have played an important role in the diversification of the clade, or speciation rates are variable in the genus. Gobies may be susceptible to adaptive radiations if they have a fast rate of molecular substitutions (e.g., Coyne and Orr, 1998; Barrier et al., 2001; Kassen, 2009), which is likely due to their extremely short generation times compared to other vertebrates (Depczynski et al., 2007). Another possible cause of increased rates of molecular substitutions is fluctuations in population size (Charlesworth, 2009; Woolfit, 2009). This may occur in gobies given that some clades have variable population sizes due to vulnerability to habitat change (Munday, 2004a,b). While increased rates of molecular substitutions may enhance adaptive radiations, ecological circumstances may have also played a significant role (e.g. the extinction of a competing clade), so further examination of the ecology of gobies should provide insight into the drivers of their radiation.

Early diversification of

While several studies on the evolution of Acropora suggest that the current diversity is mainly of Pleistocene origin (Wallace, 1999; van Oppen et al., 2001; Vollmer and Palumbi, 2002), fossil findings at high latitude Eocene assemblages indicate that a complete Pliocene diversification of the genus is untenable (Wallace and Rosen, 2006). Moreover, additional fossils suggest that the genus diversified soon after it first appeared in the fossil record in samples more than 55 My old, during the Paleocene (Carbone et al., 1993; Wallace, 2008).

Our date estimations for Acropora show wide intervals of likely diversification dates, which are likely due to significant phylogenetic uncertainty. Phylogenetic uncertainty may be caused either by incomplete lineage sorting, or introgression and hybridization, both issues pertinent to acroporid corals (Van Oppen et al., 2001). Testing these is beyond the scope of the present study, but presumably a case of introgression or hybridization would cause estimated ages to be more recent, given that these mechanisms serve as “short cuts” to evolution, blurring previous divergences (Willis et al., 2006). If the emergence of some of the species are older than the present inferences suggest, it is more likely due to poor taxon sampling than hybridization events. Furthermore, hybridization events appear to be uncommon in the wild (Vollmer and Palumbi, 2002; Márquez et al., 2002; Miller and Van Oppen, 2003; Wolstenholme, 2004), despite the potential being evident in laboratory and genetic analyses in Acropora and other scleratinian groups (Budd and Pandolfi, 2004; Van Oppen et al., 2004; Combosch et al., 2008). As such, it is worth noting that even rare hybridization events, particularly in the presence of subsequent introgression, can be evolutionary significant events (Abbott et al., 2013). Incomplete taxon sampling is also likely to mask more recent radiations. However, fishes of Gobiodon use coral species from different well-supported clades in the present results, so it is unlikely that there is a hidden pattern whereby Gobiodon fishes cospeciated with more recent species of corals that are unrepresented or unresolved in our analyses.

Although they are probably rare events, the phenomena of hybridization and introgression have resulted in surprising genetic similarity among congeners, increasing the amount of uncertainty in phylogenetic inferences and restraining these analyses from being useful to unwind the taxonomy of the genus (Van Oppen et al., 2001; Richards et al., 2008). Modern molecular dating techniques using Bayesian analysis are affected in this circumstance because they account for phylogenetic uncertainty in the estimation of dates (Drummond et al., 2006), so high genetic similarity among the samples leads to an inflation in the uncertainty incurred. This reduces the power of a dates comparison analysis to find significant differences, as the amount of overlap of the estimates of posterior distributions tends to increase. The present result, for instance, provides a broad uncertainty in the date estimates, yet shows a significant difference between the estimates for Gobiodon and its associates of Acropora.

Mitochondrial and nuclear markers provide contrasting phylogenetic relationships for Acropora (Richards et al., 2008). This could be expected due to the higher levels of introgression in mitochondrial genes than those seen in nuclear genes, as well as to the generally slow rate of evolution of corals (Van Oppen et al., 1999; Shearer et al., 2002). It is also likely that mitochondria had different evolutionary paths, causing differences in topologies and the incongruence with nuclear data on the emergence times of Acropora digitifera and A. loripes. Within these species nuclear and mitochondrial mean age estimates differ by around 2–3-fold and were significantly different given that the posterior distributions did not overlap. Although this lack of congruence impedes phylogenetic inferences of the genus, >80% of the taxa had overlapping age intervals between the two genetic datasets, precluding any further inference of difference between the mitochondrial and nuclear datasets.

A tree construction model that excludes extinction was superior to a model including extinction, which may reflect the known mechanisms of corals to reduce the propensity for extinction (Kenyon, 1997; Richards et al., 2008). These mechanisms include having genets of indefinite longevity and the capacity to hybridize, albeit rarely. The feature of avoiding extinction is congruent with a slower rate of evolution, older clade age, and a greater species diversity of Acropora compared to Gobiodon. Therefore, the evolutionary patterns of the two associates are in stark contrast; Acropora has had a slow, steady evolutionary dynamic over an extended period of time while Gobiodon has a variable, much younger evolutionary dynamic.

Conclusions

The present results indicate that the process of diversification of host corals is unlikely to have an influence on the evolution of obligate coral-dwelling Gobiodon. Instead, the present diversity of corals may have been important throughout the evolution of Gobiodon. Fishes of the genus Gobiodon comprise two major and four minor clades. The genus Gobiodon is also a substantially younger taxon than the Acropora corals they inhabit, and has a faster, more variable evolutionary history. This inference refutes the hypothesis of cospeciation. While the ecological association between Gobiodon and Acropora is likely to be a mutualism, it is asymmetrical on an evolutionary time scale. Nevertheless, the ecological attributes of the association are likely to influence the evolution of Gobiodon, which is a topic yet to be studied in detail. This is particularly relevant at a time when coral reef ecosystems are increasingly impacted by a range of anthropogenic disturbances and threats leading to enhanced risk of extinction for both Acropora corals and the diverse assemblages of fishes they support.