A world of taxonomic pain: cryptic species, inexplicable host-specificity, and host-induced morphological variation among species of Bivesicula Yamaguti, 1934 (Trematoda: Bivesiculidae) from Indo-Pacific Holocentridae, Muraenidae and Serranidae

The taxonomy of species of Bivesicula Yamaguti, 1934 is analysed for samples from holocentrid, muraenid and serranid fishes from Japan, Ningaloo Reef (Western Australia), the Great Barrier Reef (Queensland), New Caledonia and French Polynesia. Analysis of three genetic markers (cox1 mtDNA, ITS2 and 28S rDNA) identifies three strongly supported clades of species and suggests that Bivesicula as presently recognized is not monophyletic. On the basis of combined morphological, molecular and biological data, 10 species are distinguished of which five are proposed as new. Bivesicula Clade 1 comprises seven species of which three are effectively morphologically cryptic relative to each other; all seven infect serranids and four also infect holocentrids. Bivesicula Clade 2 comprises three species of which two are effectively morphologically cryptic relative to each other; all three infect serranids and one also infects a muraenid. Bivesicula Clade 3 comprises two known species from apogonids and a pomacentrid, and forms a clade with species of Paucivitellosus Coil, Reid & Kuntz, 1965 to the exclusion of other Bivesicula species. Taxonomy in this genus is made challenging by the combination of low resolving power of ribosomal markers, the existence of regional cox1 mtDNA populations, exceptional and unpredictable host-specificity and geographical distribution, and significant host-induced morphological variation.


Introduction
A key component of the understanding of biodiversity of any group is geographical distribution. In their review of the biodiversity of trematodes of fishes of the Indo-west Pacific (IWP), Cribb et al. (2016) found that 'understanding of the geographical distribution of trematode species in the IWP is especially wanting'. They pointed out that records for the region are heavily concentrated in just a few regions (especially off Hawaii, India, Japan and the Great Barrier Reef) and that 87% of known species have been reported no more than five times. The lack of records over range is exacerbated by the lack of molecular testing of parasite identity, in a group seemingly especially prone to cryptic speciation (Poulin, 2011;Pérez-Ponce de León and Poulin, 2018). Recently, a handful of studies has begun to redress this ignorance with combined morphological and molecular studies of specific taxa from widely separated IWP sites. These studies show that species from multiple families can be interpreted as geographically widespread, although often with distinct regional populations (McNamara et al., 2014;Bray et al., 2018;Bray et al., 2022;Cutmore et al., 2021;Huston et al., 2021). Less frequently, there has been evidence of separate closely related species in the same fish species across range Martin et al., 2018). These studies, especially the recognition of regional genetic variation, have crystallized problems on the basis of species recognition, leading Bray et al. (2022) to propose objective criteria for recognition of trematode species over range in the light of molecular data. This approach is applied here to trematodes of the family Bivesiculidae Yamaguti, 1934, a group for which molecular approaches have been little used (Trieu et al., 2015;Atopkin et al., 2017) and not at all over range. The Bivesiculidae is a small family of exclusively fish parasites occurring in a wide range of marine fish families (Cribb, 2002). This study aims to differentiate IWP species objectively and to characterize their host-specificity and geographical distribution.

Sample collection and morphology
Fishes were collected by spear, seine net, line or purchased from fish markets at localities off Australia, Japan, French Polynesia and New Caledonia (Table 1). Digeneans were collected from freshly killed fish as described by Cribb and Bray (2010), fixed by being pipetted into nearly boiling saline, and immediately preserved in either formalin (early work) or 80% ethanol (more recent work). Some individual worms preserved in 80% ethanol were processed for both morphological and molecular analyses (hologenophores sensu Pleijel et al., 2008) and many specimens were also designated as paragenophores (specimens collected from the same host individual as sequenced specimens).
Whole mounts for morphological analysis were stained with Mayer's haematoxylin, dehydrated in a graded ethanol series, cleared in methyl salicylate and mounted in Canada balsam. Drawings were made using a drawing tube on an Olympus BX-53 (Tokyo, Japan) compound microscope. Morphometric data were taken with the same microscope with cellSens Standard imaging software. The long axis of the testis and the ovary rarely aligned with that of the body. To capture the full size of these organs the maximum detectable length and the maximum width perpendicular to that length were measured.
PCR for the ITS2 rDNA and 28S rDNA regions was performed with a total volume of 20 μL, consisting of 5 μL of 5 × MyTaq reaction buffer (Bioline), 0.75 μL of each primer, 0.25 μL of Taq polymerase (Bioline MyTaq™ DNA polymerase), and 2 μL of DNA template, made up to 20 μL with Invitrogen™ ultraPURE™ distilled water. PCR for the cox1 mtDNA region was performed with a total volume of 20 μL, consisting of 5 μL of 5 × MyTaq reaction buffer, 2 μL of each primer, 0.25 μL of Taq polymerase and 4 μL of DNA template, made up to 20 μL with Invitrogen™ ultraPURE™ distilled water. Amplifications were carried out on a Takara TP-650 PCR thermocycler (Beijing, China). The following profile was used to amplify the ITS2 rDNA region: an initial single cycle of 95°C denaturation for 3 min, 45°C annealing for 2 min, 72°C extension for 90 s, followed by 4 cycles of 95°C denaturation for 45 s, 50°C annealing for 45 s, 72°C extension for 90 s, followed by 30 cycles of 95°C denaturation for 20 s, 52°C annealing for 20 s, 72°C extension for 90 s, with a final 72°C extension for 5 min. The following profile was used to amplify the 28S rDNA region: an initial 95°C denaturation for 4 min, followed by 30 cycles of 95°C denaturation for 1 min, 56°C annealing for 1 min, 72°C extension for 2 min, followed by a single cycle of 95°C denaturation for 1 min, 55°C annealing for 45 s, with a final 72°C extension for 4 min. The following profile was used to amplify the cox1 mtDNA region: an initial 94°C denaturation for 3 min, followed by 40 cycles of 94°C denaturation for 30 s, 50°C annealing for 30 s, 72°C extension for 30 s, with a final extension at 72°C for 10 min. Cycle sequencing of amplified DNA was carried out at the Australian Genome Research Facility with ABI Big Dye™ v.3.1 chemistry following the manufacturer's recommendations, using the same primers used for PCR amplification as well as the additional 28S primers 300F [5 ′ -CAA GTA CCG TGA GGG AAA GTT G-3 ′ ; (Littlewood et al., 2000)] and ECD2 [5 ′ -CCT TGG TCC GTG TTT CAA GAC GGG-3 ′ ; (Littlewood et al., 1997)]. Geneious® version 10.2.6 (Kearse et al., 2012) was used to assemble and edit contiguous sequences.
ITS2 rDNA and cox1 mtDNA sequence data generated during this study were each aligned with MUSCLE in MEGA 11 (Tamura et al., 2021) using UPGMA clustering for iterations 1 and 2. ITS2 rDNA sequences generated during this study were aligned with relevant sequences available on GenBank (KR092219-22, Trieu et al., 2015). The ends of each ITS2 rDNA fragment were trimmed for a final dataset of 453 base positions. The cox1 mtDNA alignment (474 base positions) was transferred to Mesquite v.3.31 (Maddison and Maddison, 2021), translated (Translation Table 9: Echinoderm/flatworm mitochondrial code) and inspected for internal stop codons. All codon positions in the cox1 mtDNA dataset were evaluated for substitution saturation, as well as non-stationarity caused by base composition bias. Substitution saturation was assessed using the 'Test of substitution saturation by Xia et al.' function (Xia et al., 2003;Xia and Lemey, 2009) as implemented in DAMBE v. 7.2 (Xia, 2018); no significant substitution saturation was detected. Non-stationarity was assessed using the χ 2 function in PAUP v. 4.0 (Swofford, 2002); significant non-stationarity was not detected. Thus, all codons in the cox1 mtDNA dataset were used in downstream analyses. Pairwise differences were estimated for each dataset using the following conditions: 'Variance Estimation Method = None', 'Model/Method = No. of differences' and 'Substitutions to Include = d: Transitions + Transversions' and 'Gaps/Missing Data Treatment = Pairwise deletion'. Unrooted Neighbour joining analyses were conducted using MEGA 11 for each dataset to explore species boundaries, with the following parameters: 'Model/Method = No. of differences', 'Substitutions to Include = d: Transitions + Transversions', 'Gaps/Missing Data Treatment = Pairwise deletion' and 'Rates among Sites = Gamma Distributed'. Nodal support was estimated by performing 10 000 bootstrap replicates.
The partial 28S rDNA sequences generated during this study were aligned with sequences of related bivesiculids from  (Miller et al., 2010), with ClustalW sequence weighting and UPGMA clustering for iterations 1 and 2. The resultant alignment was refined using MESQUITE (Maddison and Maddison, 2021) and the ends of the alignment were trimmed for a final length of 1340 base positions. Bayesian inference and maximum likelihood analyses of the 28S rDNA dataset were conducted to explore relationships among these taxa. Bayesian inference analysis was performed using MrBayes version 3.2.7a (Ronquist et al., 2012) run on the CIPRES portal and maximum likelihood analysis using MEGA 11. The best nucleotide substitution model was estimated using jModelTest version 2.1.10. The Akaike Information Criterion (AIC) and Bayesian Information Criterion (BIC) predicted the TIM3 + Γ model as the best estimator; Bayesian inference and maximum likelihood analyses were conducted using the closest approximation to this model. Nodal support in the maximum likelihood analysis was estimated by performing 1000 bootstrap pseudoreplicates. Bayesian inference analysis was run over 10 000 000 generations (ngen = 10 000 000) with two runs each containing four simultaneous Markov chain Monte Carlo (MCMC) chains (nchains = 4) and every 1000th tree saved. Bayesian inference analysis used the following parameters: 'nst = 6', 'rates = gamma', 'ngammacat = 4', and the prior parameters of the combined dataset were set to 'ratepr = variable'. Samples of substitution model parameters, and tree and branch lengths were summarized using the parameters 'sump burnin = 3000' and 'sumt burnin = 3000'. Species of the Transversotrematidae were designated as outgroup taxa (KX186733, KX186736 and KX186743, Cutmore et al., 2016).

Results
In total, 188 individuals of 15 species of Holocentridae and 807 individuals of 45 species of Serranidae from nine localities in the Indo-West Pacific were examined. In addition, bivesiculids from 275 individuals of 28 New Caledonian serranids collected by Justine et al. (2010) were examined. Infections of species of Bivesicula were found in 43 parasite/host/locality combinations. None of the species concerned were ever detected in any of thousands of individuals of other families of fishes, except for one species infecting a muraenid. Although some combinations of the specimens are morphologically distinctive, species boundaries and phylogenetic relationships were also explored with molecular data.
cox1 mtDNA analysis cox1 mtDNA data were generated, as far as possible, for at least two individuals of every available host/locality combination. A total of 74 cox1 mtDNA sequences were aligned without ambiguity or alignment gaps. An arbitrary criterion of 10 base position differences was set as potentially informative for the recognition of species or populations; these groups were treated initially as operational taxonomic units (OTUs). Differences at <10 base positions were interpreted as inconsequential intraspecific variation. On the basis of this criterion, 15 OTUs were detected (Fig. 1); the OTUs and final species level hypothesis are both labelled on the tree. Replication within these OTUs ranges from one (no replication) to 20. Three major clades (Clades 1, 2 and 3) are distinguished, two with high nodal support and one (Clade 2) with poor support; Clade 1 includes nine OTUs, Clade 2 includes four OTUs and Clade 3 includes two OTUs.

ITS2 rDNA analysis
ITS2 rDNA sequences were generated for, where possible, a minimum of two representatives from each cox1 OTU. The analysis of 41 sequences (Fig. 2) distinguished the same three major clades as the cox1 analyses, each with strong nodal support. Within Clade 1 (relating to nine cox1 OTUs), variation was exceptionally low, not exceeding four base positions. These distinctions are consistent with only four weakly differentiated OTUs. Within Clade 2, three well-supported clades were distinguished at 8-10 base positions. Within Clade 3, the two cox1 OTUs differed at a single base position.

28S rDNA analysis
28S rDNA sequences were generated for, where possible, the same two representatives of each cox1 OTU that were sequenced for ITS2 rDNA. Twenty-six sequences were analysed with other previously reported bivesiculid 28S rDNA sequence data. The results of this analysis (Fig. 3) again identified Bivesicula species as forming three highly supported clades. Sequences of species of Clades 1 and 2 differed at 27-29 base positions whereas those in Clade 3 differed from those in Clades 1 and 2 at 78-81 base positions. Clade 1 lineages showed only a low level of variation at a maximum of 4 or 5 base positions, which gave almost no useful resolution of OTUs. Clade 2 distinguished the same three OTUs as seen in the cox1 mtDNA and ITS2 rDNA analyses but with lower nodal support. Clade 3 formed a strongly supported clade

Morphology
A total of 331 slide specimens, including 30 hologenophores and 108 paragenophores were prepared and examined. These specimens corresponded to all of the OTUs referred to above except that there is no specimen relating to the single sequence from Epinephelus fasciatus (Forsskål) from off Bali. In addition to specimens corresponding to the molecular data set, multiple hostlocality combinations for which no sequences were available were examined. These specimens were measured systematically, and drawings were made of at least one individual of each of the 15 OTUs identified by cox1 mtDNA sequencing. Specimens relating to the three clades of species of Bivesicula identified by phylogenetic analyses are easily distinguished by the combination of the form of the pharynx which is large and robust in Clade 1 and proportionally far smaller in Clades 2 and 3. The two species of Clade 3 are distinguished from those of Clade 2 by their relatively small and rounded bodies with intestinal caeca exceeding the testis posteriorly.
Within Clade 1 there are four moderately distinct morphotypes with more subtle distinctions enabling the marginal recognition of a fifth. In Clade 2, there are two clearly distinct morphotypes, one of which divides imperfectly on the basis of the anterior extent of the vitelline follicles. In Clade 3 there are two clearly distinguishable morphotypes corresponding to the two previously described species, B. neglecta and B. unexpecta.

Synthesis
An integrative approach to the delineation of trematode species was adopted, considering morphology, molecular data and biology (principally host-specificity), following the species recognition criteria proposed by Bray et al. (2022). These criteria require reciprocal monophyly inferred from at least the most informative molecular marker plus either morphological distinction or biological distinction in the form of differing host distributions. Application of these criteria leads us to recognize 10 species in the collection of bivesiculids from holocentrids, muraenids and serranids. Recognition of almost all of these species can be considered difficult in that the morphological differences between them are limited. This hypothesis requires the recognition of several combinations of actually or nearly cryptic species.
Here the broad overall case for the recognition of 10 species (nine named) is made, because the argument requires reference to multiple species simultaneously.
One species, Bivesicula nana n. sp., has no molecular data but is arguably the most morphologically distinctive of those considered here on the basis of being far smaller than all the other morphologically comparable species. On the basis of its overall form, it can be predicted that it will prove to belong to Clade 1.
The nine sequenced species are all easily distinguished by cox1 data but far less so, or in some cases not at all, by ITS2 rDNA and 28S rDNA data. The two ribosomal markers are effective in distinguishing the three major clades of Bivesicula species and the three species of Clade 2 recognized, but are completely ineffectual in delineating the species of Clade 1 suggested by cox1 data. In this context, it is noteworthy that the two previously described Clade 3 species for which new sequence data are reported (B. neglecta and B. unexpecta), are morphologically clearly distinct and infect a non-overlapping range of hosts (apogonids and a pomacentrid, respectively) yet differ at only a single base position in the ITS2 rDNA alignment and three base positions in the 28S rDNA alignment. Thus, as found recently for certain Aporocotylidae , ITS2 rDNA and 28S rDNA sequences are uninformative for the distinction of some combinations of bivesiculid species. The considerations below are thus based on inferences drawn from morphology, cox1 sequences and biology (host-specificity). Reference to Figs 2 and 3 demonstrates that some of the species discussed below
(especially those of Clade 2) are distinguished by ITS2 rDNA and 28S rDNA sequences, but the taxonomic hypothesis does not depend on those data. Clade 1 here is interpreted as representing seven species if the morphologically comparable but unsequenced B. nana n. sp. is included. Of the six sequenced putative species, relating to nine OTUs, the form from E. fasciatus from Bali is represented by a single highly distinctive sequence (differing from all other taxa at a minimum of 50 base positions in the cox1 alignment); however, there is no morphological specimen available, and this form is consequently only noted here as an indication of uncharacterized richness. The five remaining species are morphologically similar and incorporate substantial difficulty in delineation, although they differ from each other at a minimum of 21 cox1 base positions. Three distinguishable morphotypes were detected, relating to that of B. claviformis, B. cephalopholicola n. sp. and B. obovata Shimazu and Machida, 1995. Of these, B. claviformis is considered actually or nearly morphologically cryptic with respect to two further species, B. sheni n. sp. and B. polynesiensis n. sp.
Samples from serranids from Okinawa are here interpreted as representative of B. claviformis, the type species for the genus. These specimens form a well-supported cox1 clade with samples from serranids from Ningaloo Reef (Western Australia) and from holocentrids and serranids from the GBR (Queensland). They differ at a maximum of 23 base positions in the cox1 alignment and relate to OTUs 4, 5, 6 and 7; these differences are here interpreted as reflecting intra-specific geographical populations. The entire clade differs from all other taxa in the analysis at a minimum of 37 base positions. The morphology of specimens corresponding to the sequenced specimens is broadly uniform, but does incorporate variation in body size, shape, the appearance of the vitelline follicles and their anterior extent. This is the only species represented in Clade 1 in which the cox1 OTUs do not correspond directly to just one of the recognized species.
Bivesicula sheni n. sp. is described as a new species explicitly cryptic relative to B. claviformis and presently known only from the GBR (at both Heron Is. and Lizard Is.). It has been collected and sequenced from two species of Sargocentron Fowler (Holocentridae) and six species of Epinephelus Bloch (Serranidae). Twenty cox1 sequences differ at only 0-4 base positions. As for B. claviformis, the morphology of specimens corresponding to the sequenced specimens is broadly uniform, but again incorporates noticeable variation in body size, shape, the appearance of the vitelline follicles and their anterior extent. There is also evidence of host-induced phenotypic variability in this species relating to infection of holocentrids and serranids. No basis for the morphological distinction of this species from B. claviformis, which occurs in an overlapping range of fishes on the GBR, was identified; indeed, the two species have been sequenced from the same individual Sargocentron caudimaculatum (Rüppell). However, in the cox1 analysis, this species differs from those identified as B. claviformis at 39-46 base positions. In addition, samples identified as B. sheni n. sp. form a wellsupported clade with samples identified as the previously described species B. obovata. Bivesicula obovata is morphologically clearly distinguishable from both B. claviformis and B. sheni n. sp. on the basis of its massive body form. The criteria for the recognition of species proposed by Bray et al. (2022) constrain us to recognize B. sheni n. sp. as a species separate from B. claviformis because, despite indistinguishable morphology and overlapping distribution and hosts, the topology of their relationships based on cox1 data does not identify them as sister taxa, but as separated by another clearly distinct species. On this basis, B. claviformis and B. sheni n. sp. can presently be distinguished where they co-occur only by sequence data (cox1). Therefore, the designated holotype for the new species is a hologenophore and the paratypes are all hologenophores and paragenophores.
Bivesicula polynesiensis n. sp. is proposed for samples from two species of Holocentridae and two of Serranidae from the Austral, Gambier and Society Archipelagos in French Polynesia. Corresponding sequences form a strongly supported cox1 clade with a minimum of 38 base positions differences relative to all other taxa. There is no population structure between specimens from the three French Polynesian collecting localities which are separated by up to 1600 km. This species closely resembles B. claviformis and B. sheni n. sp. to the point of being effectively cryptic, however it does not form a monophyletic clade with either species. Relationships between B. cephalopholicola n. sp., B. claviformis, B. sheni n. sp., B. obovata and B. polynesiensis n. sp. are poorly resolved overall, but as for the recognition of B. sheni n. sp., the close relationship between B. sheni n. sp. and B. obovata precludes recognition of B. polynesiensis n. sp. as forming a clade with B. claviformis or B. sheni n. sp. Slight tendencies to morphological distinction of B. polynesiensis n. sp. relative to B. claviformis and B. sheni n. sp. are discussed following the species account, but the three are really very similar. Complicating morphological differentiation is the fact that, as for B. sheni n. sp., there was clear evidence of host-induced phenotypic variation between samples from holocentrids and serranids. Strikingly, the variation (a tendency to a more rounded body in specimens from holocentrids than in those from serranids) is comparable to that seen for B. sheni n. sp. so that what can be interpreted as different species have convergent morphology depending on the family of hosts from which they are recovered. On present indications, this species can probably be identified reliably on the basis of its B. claviformis-like morphology in the context of infection of holocentrids and serranids in French Polynesia. However, it is clear that molecular data are necessary for identification in any other circumstances. To minimize any ambiguity associated with the description of B. polynesiensis n. sp., it is described with a hologenophore as the holotype and hologenophores and paragenophores as the paratypes.
Bivesicula cephalopholicola n. sp. is proposed for specimens from three species of Cephalopholis Bloch & Schneider from Lizard Is. This form is distinct from all other taxa at a minimum of 44 cox1 base positions. The specimens broadly resemble B. claviformis, B. sheni n. sp., B. nana n. sp. and B. polynesiensis n. sp., but they show some morphological distinction in their relatively small cirrus-sacs, biologically in their restriction to species of the genus Cephalopholis and, on the GBR, occurring only at Lizard Is., not at Heron Is. Neither B. claviformis nor B. sheni n. sp. has been detected in species of Cephalopholis at Lizard Is. although one infection of B. claviformis was found in C. argus at Okinawa. Specimens from Cephalopholis boenak (Bloch) from New Caledonia are identified as this species on the basis of their morphology and host distribution.
Clade 2 is here interpreted as comprising three species, Bivesicula gymnothoracis Shimazu and Machida, 1995, Bivesicula palauensis Shimazu and Machida, 1995 and Bivesicula novaecaledoniensis n. sp. The two previously described species are readily morphologically distinguishable by general body features. New specimens identified as B. palauensis were collected from Variola albimarginatus Baissac and Epinephelus areolatus (Forsskål) from Okinawa (with molecular data) and from Epinephelus morrhua (Valenciennes) from New Caledonia (morphology only). These samples are consistent with the original (and only previous) description of this species from the Philippines (Shimazu and Machida, 1995). Specimens identified as B. gymnothoracis were recollected from the type host, Gymnothorax kidako (Muraenidae), from Minabe, Japan and also, surprisingly, from the completely unrelated E. fasciatus (Serranidae) from Minabe and Okinawa. New specimens from G. kidako are far larger than almost all those from E. fasciatus but whether this relates to a host-induced morphological distinction, or the limited sample size is not known. Specimens from G. kidako and E. fasciatus from Minabe (OTU 13) incorporated variation at only 0-4 base positions in the cox1 alignment; samples from E. fasciatus from Okinawa (OTU 12) differed from the Minabe samples at 10-13 base positions in the cox1 dataset, a distinction interpreted as intraspecific geographical variation. Specimens from Epinephelus chlorostigma (Valenciennes) from New Caledonia are here described as a new species. This species resembles B. gymnothoracis closely, differing unreliably in only one detected morphological character, but the two species differ at 81-85 base positions in the cox1 alignment, 8-9 base positions in the ITS2 rDNA alignment, 13 in the 28S rDNA alignment, and are not sister taxa. Table 1 summarizes the geographic distribution of the 10 species on the basis of the specimens reported here. Notably, four of the 10 have been found from off Lizard Island, northern GBR. Table 2 summarize the host distribution. Notably, seven of the species have been found in E. fasciatus. All 10 species infect serranids; four also infect holocentrids and one also infects a muraenid.  (189), or 13.4-21.3 (16.7) % body length, from anterior extremity, usually almost at level of anterior arms of excretory vesicle but sometimes leaving tips clearly exposed, to 219-508 (370), or 25.4-36.5 (32.2) % body length, from posterior end of body, usually terminating distinctly anterior to ends of intestinal caeca but sometimes extending almost to their termination; total vitelline field 367-815 (584) long or 45.2-56.3 (51.2) % of body length. Uterus passes to close to posterior extremity before passing anteriorly to open at common genital pore. Eggs 68-90 × 35-49 (44 × 77). Excretory vesicle V-shaped; arms pass latero-dorsally to testis, terminating distinctly posterior to posterior margin of pharynx, 124-205 (161), or 11.5-18.3 (14.4) % body length, from anterior extremity. Excretory pore terminal.

Remarks
This species is the type species for the genus and, as is commonly the case, it is by far the most frequently reported species in the genus. There remain some problems in the understanding of the status of this species. The type-host was reported as Seriola quinqueradiata Temminck & Schlegel (Carangidae) from the Inland Sea of Japan (Yamaguti, 1934). Apparently, no carangid has been reported as a host for a bivesiculid since. Thus, if the species is actually one that shows high specificity for S. quinqueradiata, then it may be that the species has not been re-reported since its first description; however, that seems unlikely and probably the original infection was exceptional. Yamaguti (1938) described a second species, B. epinepheli Yamaguti, 1938, from Epinephelus akaara (Temminck & Schlegel), also from the Inland Sea. Although Yamaguti (1938) stated that this species differs notably from B. claviformis in the position of the pharynx, the differences appear inconsequential and the two species were synonymized by Cribb et al. (1994), an action accepted by Shimazu (1994) and Shimazu and Machida (1995). Bivesicula xishaensis Gu and Shen, 1983 was described on the basis of a single specimen from Epinephelus fasciatus collected in the Xisha Islands, off China (Gu and Shen, 1983). The description lacks some important details (e.g. the anterior extent of the excretory vesicle) but generally resembles B. claviformis and it was synonymized with that species by Cribb et al. (1994).
The specimens (hologenophores and paragenophores) corresponding to sequences generated in this study from C. argus, E. fasciatus and E. merra from Okinawa are broadly consistent with the original descriptions of B. claviformis (and B. epinepheli) from the same area. More recent records of B. claviformis by Shimazu and Machida (1995) and Kuramochi (2018) from Japanese waters have been overwhelmingly from serranids so that, despite the issue with the identity of the type host for B. claviformis, the molecular identity of this species can probably now be considered established with reasonable confidence. The cox1 sequences from Okinawan serranids form a strongly supported clade with sequences from samples from E. fasciatus from Ningaloo Reef and from E. merra, S. caudimaculatum and S. spiniferum from Lizard Is. Variation between geographic populations of this clade are at 17-23 base positions in the cox1 dataset, a level consistent with what has been interpreted as intraspecific geographic variation for other trematodes such as species of Hurleytrematoides (see McNamara et al., 2014) and Preptetos (see Bray et al., 2022). Differences between these localities in the ITS2 and 28S rDNA alignments are no greater than 1 and 2 base positions, respectively.
The specimens of this species are by no means identical; there is variation in body shape, size of the vitelline follicles and appearance of the pharynx. However, no regional pattern of variation was detected (variation in specimens from a single locality was as great as that between localities) and the specimens are broadly consistent with each other. All these forms are thus interpreted as B. claviformis.
The evidence discussed above suggests that B. claviformis occurs widely in the Indo-Pacific. However, the molecular analyses also show the presence of two other strongly supported independent cox1 clades from French Polynesia and the GBR, respectively, for which the corresponding specimens have highly similar morphology. Given the topology of the relationships between the various clades of species of Bivesicula, those clades are recognized as representing distinct species. The description of those species calls into question the reports of B. claviformis as being widespread across the IWP. Reports from China (Gu and Shen, 1983), the Red Sea (Nagaty, 1948), Fiji (Manter, 1961;Rigby et al., 1997;Nahhas et al., 2004) and Indonesia (Fischthal and Kuntz, 1965) all require genetic characterization and further morphological study to explore their status. It is unlikely to be a coincidence that the area most heavily sequenced for specimens resembling B. claviformis, the GBR, is that for which there is evidence of sympatric cryptic species. It is thus possible, perhaps probable, that sympatric cryptic species will also be found elsewhere.

Remarks
Because of the difficulty of distinguishing this species from B. claviformis (see below), the species is based entirely on type specimens that are hologenophores and paragenophores. Measurements given above all relate to paragenophores. However, in the prevalence listing above, unsequenced specimens have been tentatively identified as B. sheni n. sp. to give an indication of the overall prevalence detected. The specimens are tentatively listed as B. sheni n. sp. on the basis that overall sequencing shows that species to be far more common in GBR fishes than is B. claviformis, which has only been found so far in the northern GBR. Notably, three of the 13 host/locality combinations (E. cyanopodus and S. rubrum from Heron Is. and E. ongus from Lizard Is.) have never been sequenced. Supplementary Table 2 lists measurements of 96 adult whole-mount specimens that are tentatively assigned to this species; the 96 includes 30 paragenophore specimens. Figure 6A shows the distribution of body length relative to width for specimens interpreted as B. sheni n. sp. distinguished by host family (Holocentridae or Serranidae); although there is overlap in body shape of specimens from the two fish families, there is a clear tendency for specimens from holocentrids to be broader than those from serranids. In addition, 67 of 80 specimens measured from serranids are longer than the longest of 16 specimens from holocentrids suggesting strongly that total size is also affected by host identity. It was thus a surprise when sequence data indicated that samples from GBR holocentrids and serranids were the same species; their appearance at the time of collection together with preconceptions with respect to host-specificity led to an expectation that they represented different species.
In the combination of the intestinal caeca not extending posterior to the testis, the vitelline follicles not reaching the posterior margin of the pharynx or extending significantly anterior to the excretory vesicle, and the pharynx being larger than the eggs, this species is clearly distinct from all existing species of Bivesicula except for B. claviformis and B. lutiani. Of these, Bivesicula lutiani is a problematic species. It was described on the basis of 23 specimens (a robust sample) from one of six Lutjanus kasmira examined from off the Xisha Islands, Japan by Gu and Shen (1983). The host record is intriguing because species of Lutjanus, including L. kasmira, have been heavily studied for their digenean parasites without other reports of species of Bivesicula. Thus, the infection is either atypical (a species normally infecting another fish species) or exceptional (genuine, but isolated). Gu and Shen (1983) observed that the species closely resembles B. claviformis Yamaguti, 1934 andB. epinepheli Yamaguti, 1934 (now a synonym of B. claviformis). It was distinguished from the former by the size of body, the width and position of the cirrus sac, and the position of the testes, and from the latter in the shape and the size of the body, the position of ovary and the size of the eggs. In light of the level of variation seen in other species of Bivesicula, it seems likely that none of these characters reliably delineates B. lutiani. The status of this species requires further attention.
The key comparison for B. sheni n. sp. is relative to the type species, B. claviformis. Samples here identified as B. sheni n. sp. differ from those of B. claviformis (reported above) at 39-45 base positions in the cox1 dataset. Both taxa form strongly supported clades. Most significantly, in terms of the recognition of B. sheni n. sp. as a new species, is the observation that its closest relative in analysis of cox1 sequences is B. obovata (characterized below), from which it differs clearly in morphology and with which it forms a strongly supported clade. Bivesicula sheni n. sp. and B. obovata infect an overlapping range of holocentrid and serranid fishes. When ITS2 and 28S rDNA sequences are considered, there are almost no differences between B. claviformis, B. sheni n. sp. and B. obovata. Corresponding ITS2 rDNA sequences may be identical between the three and never differ by more than 2 base positions. Corresponding 28S rDNA sequences differ between B. claviformis and the other two forms at 2-4 base positions and between B. sheni n. sp. and B. obovata at 0-2 base positions. None of the differences for either marker results in the recognition of clades distinguishing any of these forms which occur sympatrically on the GBR. In short, B. claviformis and B. sheni n. sp. are morphologically indistinguishable, a problem exacerbated by the considerable intraspecific variation in both, including host-induced morphological variation (dependent on host family) for B. sheni n. sp. Charting of numerous morphometric characters for the two forms failed to detect any reliable basis for distinction. The vitelline follicles of B. claviformis tend to be slightly more anteriorly extensive than those of B. sheni n. sp., but the distinction is not reliable. Bivesicula sheni n. sp. is thus here recognized explicitly as a species cryptic relative to B. claviformis and is distinguished on the basis of consistent genetic differences (principally in cox1) and the fact that it does not form a clade with B. claviformis.

Remarks
Bivesicula obovata was described by Shimazu and Machida (1995) on the basis of three specimens, two mature and one immature, from a holocentrid, S. rubrum, from Okinawa, Japan. The new specimens, including some from S. rubrum, are consistent with their description, most noticeably in the distinctively broad body shape. The original description and figure reports the vitelline follicles as extending just anterior to the arms of the excretory vesicle as opposed to just short of them in the new specimens, but this distinction is considered insignificant. The specimens of this species were noted as being probably distinct from all others at the time of collection on the basis of their massive body shape. The species has been found by us only at Heron Is. where it has never been common in any of the four hosts in which it has been detected. Only two intact paragenophore specimens and two hologenophores (one front end, one back end) were available to consider the morphology of this species; the single specimens from E. fasciatus and E. undulatostriatus were consumed entirely for molecular analysis. In terms of cox1 sequence data, this species is closest, by far, to B. claviformis with which it forms a strongly supported clade. However, Fig. 6B shows the clear distinction in body shape for the two intact GBR specimens together with that from the original description relative to 95 samples of B. sheni n. sp., including a single (non-paragenophore) specimen from S. rubrum from Heron Is. Surprisingly, despite the ease with which these forms are distinguishable on the basis of body shape, no other morphometric character was found that clearly distinguishes the two species.

Remarks
As for B. sheni n. sp., there is a clear tendency for specimens of this species from holocentrids to be less elongate than those from serranids (Fig. 6C), although there is overlap in body shape for the smaller specimens. If specimens from holocentrids and serranids for B. sheni n. sp. and B. polynesiensis n. sp. are plotted together, those from holocentrids of the two putative species are more similar to each other in body shape than to those of the same species (Fig. 6D). In addition, and again similar to the findings for B. sheni n. sp., 3 of 13 specimens measured from serranids are longer than the longest of 20 specimens from holocentrids suggesting that total size is also affected by host identity.
This species is clearly genetically close to B. claviformis, B. cephalopholicola n. sp., B. obovata, and B. sheni n. sp. In the cox1 analysis it forms a clade weakly sister to the highly supported clade of B. sheni n. sp. + B. obovata and B. claviformis. The hosts of all four putative species broadly overlap and all have been found in E. fasciatus. In terms of morphology, B. polynesiensis n. sp. is easily distinguished from B. obovata on the basis of the broad body shape of the latter. It is far more similar to B. claviformis and B. sheni n. sp. Plotting of pharynx length relative to body length suggests a tendency for the pharynx to be larger in B. polynesiensis n. sp. than in the other two species; the effect is clearest in larger specimens but there is significant overlap, so this does not create a reliable basis for distinction of the species. Thus, B. polynesiensis n. sp. is recognized as a species explicitly cryptic relative to B. claviformis and B. sheni n. sp. and designated type specimens are all hologenophores and paragenophores.

Remarks
This species has been detected in three species of Cephalopholis at Lizard Is. on the northern GBR and in one of the same species from New Caledonia. The New Caledonian identification is not supported by molecular data, but the shared host and comparable morphology lead us to identify them as the same species with confidence. Despite the examination of 64 individuals of two susceptible Cephalopholis species at Heron Is., this species has so far not been detected there, although the most commonly infected host species, C. microprion, has not been examined there. Single infections of species of Bivesicula in species of Cephalopholis have been detected in French Polynesia and at Okinawa, but, on the basis of sequence data, neither relates to this species.
The morphology of this form is clearly broadly comparable to that of B. claviformis, B. obovata, B. sheni n. sp. and B. polynesiensis n. sp. Relative to these species it is moderately morphologically, biologically and genetically distinctive. cox1 sequences of B. cephalopholicola n. sp. differ from all other species at a minimum of 45 base positions. Within Clade 1, based on cox1 and 28S rDNA, it is sister to all other species except for the form from Bali. Morphologically it is far narrower than B. obovata. It is similar in shape to B. claviformis, B. sheni n. sp. and B. polynesiensis n. sp. but has a generally proportionally smaller cirrus-sac (although there is considerable overlap in measurements) (Fig. 6E). The data for the GBR suggest that the species is restricted to the serranid genus Cephalopholis and, that within the known range, no other bivesiculids infect the same fishes. These conclusions are supported by the findings from New Caledonian serranids. Specimens identified as B. cephalopholicola n. sp. were found only in C. boenak; absence of infections of morphologically comparable specimens in large numbers of multiple species of Epinephelus supports the inferred restriction of B. cephalopholicola n. sp. to species of Cephalopholis. Bivesicula nana n. sp. (Fig. 9C) Type host Epinephelus maculatus (Bloch) Serranidae.

Remarks
This is the only species dealt with here for which molecular data are lacking. The species was found only once, in E. maculatus at Lizard Is. However, the morphology of this form is so distinctive that a new species can be proposed for it with confidence. The species generally resembles B. claviformis and its relatives and, like them, infects a serranid. The distinctiveness of the species relates almost entirely to its size. The five good quality gravid specimens range in length only from 502 to 647 μm. This small size is not unprecedented in the family. Bivesicula neglecta is reported at only 437-675 μm, B. unexpecta at only 506-633 μm and B. caribbensis Cable & Nahhas, 1962 at 654-1020 μm. However, these three species are all immediately distinct from the present form in having bodies that are nearly round (length/width ratio of <1.50 in contrast to that of 2.48-2.88 for B. nana n. sp.). None of the species with body shapes and proportions generally comparable to the present form have ever been reported in the size range reported here. Confidence in the value of this morphological distinction is strengthened by comparison with 119 specimens interpreted as representing B. sheni n. sp. Of these, 17 were immature but at least as large as the smallest gravid B. nana n. sp. and 14 were larger than all of the B. nana n. sp. The 90 gravid B. sheni n. sp. all exceeded 800 μm. These distinctions, illustrated in Fig. 6F, included specimens of both B. sheni n. sp. and B. nana n. sp. from the same individual E. maculatus so that there is no basis for host-induced morphological distinction in these differences. In addition to the difference in body size, the eggs of B. nana n. sp. are also noticeably smaller than those of any of the other species considered here. Although the eggs of B. australis, B. megalopis Shen, 1982 andB. tarponis Sogandares-Bernal &Hutton, 1959 are all reported at about the size of the eggs of B. nana n. sp., none has as small a maximum length ( just 59 μm). Again, none of these species shows any general resemblance to B. nana n. sp.

Remarks
This species was described effectively by Shimazu and Machida (1995) as a distinctive species. As such, no new description is provided, but new measurements and a figure for a new host/locality combination are provided. Bivesicula palauensis is among the largest species of Bivesicula, reported to over 3 mm, and also distinctive in its small pharynx relative to overall body size and especially in having the vitelline follicles extend to <15% of the body length from the posterior extremity. Most species in the genus have the vitelline follicles at least 20% of the body length from the posterior extremity. The only other species comparable in this character is B. neglecta, a much smaller and highly squat species that is otherwise not confusable with B. palauensis. Of six specimens of B. palauensis collected from Okinawa, two were sequenced as hologenophores, two damaged specimens were consumed entirely for sequencing, and two moderately wellfixed specimens were stained and mounted intact. From New Caledonia a single gravid specimen and two immatures from one individual E. morrhua were reported by Justine et al.
(2010) as an unidentified species of Bivesicula. All these specimens are consistent with the original description of the species, Parasitology especially in the close approach of the vitelline follicles to the posterior extremity. New specimens were found in species of Variola and Epinephelus, just as originally reported by Shimazu and Machida (1995), although from different species. Although the specimens are smaller than reported by Shimazu and Machida (1995), they are still larger than almost all the other specimens of other species of Bivesicula reported here. The new samples are confidently identified as B. palauensis, extending the known distribution from Palau to Okinawa and New Caledonia. Sequence data suggest that this species is unambiguously different from all the other species of Bivesicula for which comparable data are now available.

Remarks
This species was described effectively by Shimazu and Machida (1995) on the basis of three specimens taken from G. kidako from Fukaura, Ehime Prefecture, Japan. New specimens were obtained from G. kidako from the Minabe fish market, a little under 300 km to the north-east of Fukaura. Although the original specimens are all larger than those reported from G. kidako here (3.08-4.56 mm long v. a max. of 2.99 mm), the morphological agreement and the identification is clear. Specimens of this species were also collected from the serranid E. fasciatus (Fig. 10C and D) at both Minabe and in Okinawa. This identification came as a great surprise, given the host distinction, and was revealed initially only by the generation of identical cox1 sequences for specimens from the two fish species. The size range of the samples from the two host species is almost nonoverlapping. Nine of 10 good specimens from G. kidako are 2.42-2.99 mm long whereas 15 from E. fasciatus are 1.02-1.36 mm long. However, a single smaller specimen, 1.29 mm long detected from G. kidako, agrees well with those from E. fasciatus. Although morphologically comparable, B. palauensis and B. gymnothoracis are immediately distinguishable by the more posteriorly extensive vitellarium of the former species.

Remarks
This species falls clearly in Clade 2, with B. palauensis and B. gymnothoracis, in all molecular analyses. It is known from E. chlorostigma from New Caledonia on the basis of 35 specimens and multiple sequences and just three specimens from E. fasciatus, not sequenced, from the same locality. Despite the limited samples from E. fasciatus, their morphology is strongly consistent with that of the specimens from E. chlorostigma. This species is easily distinguishable from B. palauensis on the basis of its significantly less posteriorly extensive vitelline follicle distribution. However, it is strikingly similar to B. gymnothoracis which is surprising given the substantial molecular distinctions between the two (81-84 cox1, 7 ITS2 rDNA and 11 28S rDNA base positions differences between the two species), equivalent to the differences between multiple combinations of other clearly morphologically distinct species. Just one character showed some distinction between the twothe anterior extent of the vitelline follicles relative to the pharynx. Although measurements for this character overlap extensively for the smallest specimens, they diverge substantially for larger specimens. No basis for distinction was found in body shape, the posterior extent of the vitelline follicles, size of the pharynx, or egg size. Significantly, the two species have an overlapping host range (both infect E. fasciatus).
In the light of the criteria for the recognition of species, the failure of B. gymnothoracis and B. novaecaledoniensis n. sp. to form an exclusive clade (B. novaecaledoniensis n. sp. is consistently identified as sister to B. palauensis) means that a new species is required for these specimens. Given that they are morphologically cryptic, the designated type specimens are all paragenophores.
Justine et al. (2010) examined 61 individual E. fasciatus from New Caledonian waters for internal parasites and, for bivesiculids, found only the three specimens here identified as B. novaecaledoniensis n. sp. in three individual fishes. They also examined 38 E. maculatus and 18 E. merra which harbour B. sheni n. sp. on the GBR. These data suggest the possible absence of B. sheni n. sp., B. claviformis and B. obovata in New Caledonia, but their replacement with B. novaecaledoniensis n. sp. which has not been found on the GBR. This possible absence is surprising given B. claviformis has been reported from Fiji by Manter (1961).

Species recognition
If a simple morphological species concept was employed, it would be reasonable to recognize six species of Bivesicula in this collection from holocentrids, muraenids and serranids: B. cephalopholicola n. sp., B. claviformis, B. gymnothoracis, B. nana n. sp., B. obovata and B. palauensis. This classification would, however, have had difficulty recognizing host-induced morphological variation for some of the species. Molecular data enable the recognition of three additional largely cryptic species, B. novaecaledoniensis n. sp., B. polynesiensis n. sp. and B. sheni n. sp., suggest the presence of a fourth on the basis of a sequence of a worm from Bali, confirm the distinction of all the others except for B. nana n. sp. which remains unsequenced, and allow the recognition of host-induced morphological variation.
The final taxonomic hypothesis is heavily dependent on the interpretation of cox1 sequence data. Most importantly, combinations of samples ultimately interpreted as B. claviformis, B. sheni n. sp. and B. obovata differ at 44-47, 38-43 and 21-24 base positions in the cox1 alignment. However, they have partly identical ITS2 rDNA sequences that never differ at more than a single base position and for 28S rDNA sequences the same three combinations differ at only 1-2 base positions. On the basis of evidence relating to two previously described species it appears that ITS2 and 28S rDNA are poor markers for closely related species of Bivesicula. Bivesicula neglecta and B. unexpecta are clearly morphologically distinct and infect fishes of separate families (Pomacentridae and Apogonidae) in sympatry. They differ at just one base position in the ITS2 rDNA alignment (Trieu et al., 2015) and three base positions of 28S rDNA (this study) but new data here shows that they differ at 35 base positions in the cox1 alignment. The case of the delineation of B. sheni n. sp. and B. obovata in this study is similarly informative. The two are readily morphologically distinguishable, although in this case they infect an overlapping range of fishes in sympatry. The two forms have identical ITS2 rDNA sequences and differ by just one base position in 28S rDNA data but are clearly distinguished at 21-24 base positions in the cox1 dataset. Despite its extensive successful use to delineate species of many trematode families (see Blasco-Costa et al., 2016), it appears that ITS2 and 28S rDNA sequences are frequently inadequate markers for species delineation in this family. A comparable problem was identified recently for certain Aporocotylidae .
This study applied the criteria proposed by Bray et al. (2022) for recognition of species (reciprocal monophyly plus distinction in either morphology or host-specificity). Testing for reciprocal monophyly using only cox1 mtDNA sequences has the problem that this marker does not recombine easily and is highly prone to rapid divergence so that what might best be ultimately interpreted as regional populations may be strongly reciprocally monophyletic. It is for that reason that, for an effective basis for species recognition, Bray et al. (2022) called for the additional fulfilment of at least one of the two additional criteria. By this approach they recognized multiple monophyletic, regionally distributed populations as relating to a single lepocreadiid species (Preptetos laguncula Bray & Cribb, 1996) because, ultimately, they formed a single clade, were morphologically indistinguishable, and infected exactly the same host species. In the present system the two most problematic combinations of species are resolved on the following bases: B. sheni n. sp. and B. obovata form a single clade, infect an overlapping range of fishes but are also reciprocally monophyletic and morphologically distinguishable; B. claviformis, B. sheni n. sp. and B. polynesiensis n. sp. are morphologically indistinguishable (B. polynesiensis n. sp. is weakly distinguishable) and infect the same hosts, but each is reciprocally monophyletic with respect to each other and no combination of the three forms a monophyletic clade. Critical to this topological issue is the presence of the morphologically distinguishable B. obovata forming a clade with B. sheni n. sp.; without the presence of this species there would be a case to recognize B. claviformis, B. sheni n. sp. and B. polynesiensis n. sp. as a single widespread species. Were B. claviformis and B. sheni n. sp. to be considered a single species then the problem of cryptic species in the same sympatric hosts (on the GBR) would be removed. However, another discrepancy would arise in that B. claviformis as recognized here is a widespread species with a commensurate population structuredistinct cox1 populations in Japan, Ningaloo and the GBR distinguished at a level (17-23 bases positions) readily interpreted as intraspecific variation. The presence of one of these populations on the GBR together with a second form (here recognized as B. sheni n. sp.) that differs at 44-46 base positions in the cox1 dataset seems inconsistent with the radiation of a single species. Notably, Bray et al. (2022) reported that two species of Preptetos are each represented by two distinct cox1 populations on the northern GBR. The final distinction between the recognition of those samples of Preptetos as representing single species and the specific distinction between B. claviformis and B. sheni n. sp. recognized here is dependent on the finding that the Preptetos lineages ultimately formed monophyletic assemblages whereas B. claviformis and B. sheni n. sp. do not.
The subtleties of species recognition discussed here led to the conclusion that species of Bivesicula are so morphologically conservative that their morphology often reflects species level differences only weakly, if at all. This point is emphasized by the finding that species of Bivesicula Clade 3 form a clade with the highly distinctive genus Paucivitellosus rather than with the other clades of Bivesicula. This topology suggests that Bivesicula will ultimately need division but that the morphological distinctions between the resulting genera will be subtle. Two general explanations suggest themselves for the apparent morphological conservatism of species of Bivesicula. The first combines the possibilities that species of Bivesicula are unusually morphologically constrained and unvarying for trematodes, that this constraint is exacerbated by the absence of key digenean characters (the suckers) typically used heavily in taxonomy, and the possibility that the taxa concerned are a recently and rapidly diverging group of species. The alternative possibility is that there is little special about the group (apart from the lack of suckers) but that instead it has been studied over the range at a depth which is not yet common for trematodes of fishes of the Indo-Pacific. It is clear that the complexities of systems such as these make the issue of species recognition difficult. It is to be hoped that the stable criterion-based and system-specific approach to recognition of species used here at least has the capacity to make taxonomic proposals consistent.

Host specificity
The evolving understanding of the host specificity of bivesiculids is proving exceptionally interesting. The dominant paradigm of fish trematode specificity, at least for those of tropical marine fishes, is of oio-or stenoxenicity (Miller et al., 2011), the infection of single or closely related species. In this context it was astonishing to find clear evidence of euryxenicity (infection of phylogenetically unrelated species) for four species shared between holocentrids and serranids and one shared by a serranid and a muraenid; in each case the two families belong to different orders of fishes. The five other species in the study have been found only in serranids, but in each case sampling is probably inadequate to allow confidence that the species is restricted to just one family of fishes. The basis of this exceptional pattern of host-specificity is far from clear but must be partly based on the mode of transmission. Cribb et al. (1998) provided evidence that transmission of B. claviformis on the GBR (probably actually the new species B. sheni n. sp.) is via small fish which had ingested cercariae and in which the juvenile parasite awaited transmission in the gut. Such a transmission pattern, dependent on ingestion of small fish, is broadly consistent with the diets of a wide range of holocentrids, muraenids and serranids. However, there is no clear explanation for the general absence of infections in the many other families of piscivorous fishes. For example, on the GBR B. sheni n. sp. infects two orders of fishes, seven species of Epinephelus (Serranidae) and two species of Sargocentron
(Holocentridae). It has not been found in other well-sampled serranids such as species of Plectropomus Oken, which are voracious piscivores as shown by their rich fauna of bucephalids (Bott et al., 2013) or in species of Cephalopholis which harbour a different species of Bivesicula sympatrically (this study). These discrepancies have no clear explanation. The anomalous sharing of unrelated host families might be considered a species-specific character and thus argue against the interpretation that four species share infection of holocentrids and serranids. However, there is no reason per se to argue against a modest radiation of a lineage of Bivesicula with this exceptional host specificity. The finding that B. gymnothoracis, belonging to a clearly distinct lineage of Bivesicula species, infects both muraenids and serranids (again fishes of separate orders), supports that idea. The only comparably disjunct host distribution of GBR fish trematodes is that of Transversotrema borboleta Hunter and Cribb, 2012, infecting both chaetodontid and lutjanid fishes which are distantly related, but none of many other seemingly suitable hosts (Hunter and Cribb, 2012). The second aspect of the host-specificity of the group of Bivesicula species considered here is the remarkable role of Epinephelus fasciatus; seven of the 10 species distinguished here have been found in this fish. No aspect of the biology of E. fasciatus suggests why it should be especially susceptible to infection by bivesiculids. Strikingly, E. fasciatus is infected by two lineages of species (Clade 1 -B. claviformis, B. sheni n. sp., B. obovata, B. polynesiensis n. sp. and the form known only from Bali; Clade 2 -B. gymnothoracis and B. novaecaledoniensis n. sp.). An issue not yet much considered in the structuring of the trematode fauna of Indo-Pacific fishes is that of how widespread fish species are exploitedis it by the same or different species at different localities? There is a developing body of molecular evidence supporting widespread distributions for trematode species (e.g. Aiken et al., 2007;Bray et al., 2022;Huston et al., 2021), although increasingly populations may be divided geographically (as here for B. claviformis) (Bray et al., 2018(Bray et al., , 2022Cutmore et al., 2021). In the present study there is evidence for another pattern, speciation independent of that of the definitive hosts. This system is rendered complex by the conclusion that at least one species (B. claviformis) is widespread and shows intraspecific geographic structuring. Certainly, E. fasciatus is an especially interesting subject for further analysis over its range. In this context it is noteworthy that Kuriiwa et al. (2014) demonstrated that E. fasciatus itself has a complex structure of three cryptic lineages that may co-occur in Japanese waters. The relationship that this structure bears to the nature of their trematodes remains to be explored.

Host-induced phenotypic plasticity
The final aspect of this system that is exceptional and problematic, is the evidence that the morphology of some species of Bivesicula varies noticeably with the host in which they develop. The broad experience of trematodes of the gastrointestinal tract of tropical marine fishes, is that such phenotypic plasticity is not common; however, infection of such unrelated hosts is also not common. Perusal of the literature also suggests that the phenomenon is not reported especially frequently. There is a potential interpretative trap in this area; in the absence of molecular or experimental data, phenotypic plasticity might well go unrecognized as such and be interpreted as inter-specific variation. Regardless, there are reports relating to the phenomenon for at least 10 trematode families: Diplostomidae (Pérez-Ponce de León, 1995); Echinostomatidae (Hildebrand et al., 2015); Gorgoderidae (Bakke, 1988;Cutmore et al., 2010); Haematoloechidae (Kennedy, 1980); Haploporidae (Tonatiuh Gonzalez-Garcia et al., 2021); Heterophyidae (Elshazly et al., 2007;Fraija-Fernandez et al., 2015;Presswell and Bennett, 2020); Microphallidae (Martorelli and Ivanov, 1996); Plagiorchiidae (Blankespoor, 1974); Schistosomatidae (Machadosilva et al., 1994); Zoogonidae (Mouahid et al., 2012).
The host-induced morphological variation that is seen here has two significant facets. First, the variation occurs in the context of morphology of species of a genus in which distinguishing characters are relatively few; for a start, dimensions and relativities of oral and ventral suckers are not available for analysis. In addition to the host-induced morphological variation, there was considerable intraspecific variation not obviously associated with any biological or geographical variable. Thus, there are multiple sources of intraspecific morphological variation among these taxa, increasing the difficulty of finding reliable characters for interspecific distinction. Second, the remarkable disjunct host distributions of many of these species can easily create a false expectation of the likely presence of separate species which is reinforced by subtle distinctions in morphology. The combination of these effects, together with the complexities of the distribution of the species involved, mean that molecular data are essential for progress with this field.
Supplementary material. The supplementary material for this article can be found at https://doi.org/10.1017/S0031182022000282.
Data. The data that support the findings of this study are available from the corresponding author upon reasonable request.