Diversity, host specificity and biogeography in the Cladocorynidae (Hydrozoa, Capitata), with description of a new genus

The hydrozoan family Cladocorynidae inhabits tropical to temperate waters and comprises the two genera Pteroclava and Cladocoryne. Pteroclava lives in association with some octocorals and hydrozoans, whereas Cladocoryne is more generalist in terms of substrate choice. This work provides a thorough morpho‐molecular reassessment of the Cladocorynidae by presenting the first well‐supported phylogeny of the family based on the analyses of three mitochondrial and four nuclear markers. Notably, the two nominal genera were confirmed to be monophyletic and both morphological and genetic data led to the formal description of a new genus exclusively associated with octocorals, Pseudozanclea gen. nov. Maggioni & Montano. Accordingly, the diagnosis of the family was updated. The ancestral state reconstruction of selected characters revealed that the symbiosis with octocorals likely appeared in the most recent common ancestor of Pteroclava and Pseudozanclea. Additionally, the presence of euryteles aggregation in the polyp stage and the exumbrellar nematocyst pouches with euryteles represent synapomorphies of all cladocorynid taxa and probably emerged in their most recent common ancestor. The analysis of several Pteroclava krempfi colonies from Indo‐Pacific and Caribbean localities associated with several host octocorals revealed a high intra‐specific genetic variability. Single‐ and multi‐locus species delimitations resulted in three to five species hypotheses, but the statistical analysis of morphometric data showed only limited distinction among the clades of P. krempfi. However, P. krempfi clades showed differences in both host specificity, mostly at the octocoral family level, and geographic distribution, with one clade found exclusively in the Caribbean Sea and the others found in the Indo‐Pacific.


Introduction
Coral reefs are well known to host a great diversity of symbiotic relations, with scleractinian corals being associated with numerous other organisms, spanning all domains of life (Gates and Ainsworth, 2011;Stella et al., 2011;Hoeksema et al., 2012. Coral associates may interact with each other and the host, resulting in positive outcomes for the whole coral symbiome (Ainsworth et al., 2020), but in many cases these interactions are still poorly characterized, in particular for their taxonomic composition and ecological relevance (Gates and Ainsworth, 2011). Other benthic coral-reef invertebrates also provide habitats for a plethora of organisms (e.g., Antokhina and Britayev, 2012;Neo et al., 2015;Sch€ onberg et al., 2015;Garc ıa-Hern andez et al., 2019), even if these associations are rarely well characterized (Montano, 2020).
Octocorals represent a prominent group among nonscleractinian reef dwellers (Fabricius and Alderslade, 2001), and, due to their topological complexity and abundance, they host several associated organisms (Goh et al., 1999;Hoeksema et al., 2015). For instance, Maggioni et al., (2020a) recently reported the octocoral Bebryce cf. grandicalyx to recurrently host the hydroid Zanclea timida Puce, Di Camillo and Bavestrello, 2008 together with a suberitid sponge in the Maldives, and many other invertebrates were observed to dwell on the octocoral or sponge surface. This octocoral-hydrozoan symbiosis is not an exception, as many hydrozoan and octocoral species are reported to live in intimate associations (Bo et al., 2011). Some of these hydroids appear to be obligate symbionts of octocorals and have never been reported growing on other substrates (Puce et al., 2008a, b).
The aforementioned Z. timida, and the cladocorynid Pteroclava krempfi (Billard, 1919) are two examples of obligate alcyonancean-associated hydrozoan species (Puce et al., 2008b). Both species are known to grow on a variety of hosts: Z. timida is associated with the alcyonacean genera Paratelesto Utinomi, 1958 and Bebryce Philippi, 1841 (Puce et al., 2008b;Maggioni et al., 2020a), whereas P. krempfi associates with at least nine alcyonacean genera (Billard, 1919;Varela and Cabrales Caballero, 2010;Seveso et al., 2016Seveso et al., , 2020Maggioni et al., 2016;Montano et al., 2017). Genetic data of P. krempfi have revealed the presence of a complex of multiple cryptic species with a possible host specificity at the octocoral genus or family level Montano et al., 2017). Specifically, a Maldivian clade was associated with Paraplexaura K€ ukenthal, 1909(family Plexauridae Gray, 1859, a second Maldivian clade was found exclusively on the Alcyoniidae Lamouroux, 1812 genera Sinularia May, 1898, Sarcophyton Lesson, 1834, and Lobophytum Marenzeller, 1886, whereas a third Caribbean clade was associated with Antillogorgia Bayer, 1951(family Gorgoniidae Lamouroux, 1812. Despite the host preferences and genetic divergence, no morphological differences were detected among the three P. krempfi clades. This, along with (i) the lack of suitable genetic material of P. krempfi from the type locality (in Vietnam), and (ii) the type's host Cladiella krempfi (Hickson, 1919) and (iii) the scant information on the congeneric and hydrozoan-associated species Pteroclava crassa (Pictet, 1893), did not allow the authors to name the species.
Taxonomic confusion is also evident in the other cladocoryinid genus, Cladocoryne Rotch, 1871. All Cladocoryne species are characterized by the unique aboral branched tentacles (Bouillon et al., 2006), and the number and arrangement of aboral tentacles are key distinguishing characters for the species of this genus (Schuchert, 2003). This taxon currently contains five species, three of which [Cladocoryne travancorensis (Mammen, 1963), Cladocoryne littoralis (Mammen, 1963), Cladocoryne minuta Watson, 2005] have never been observed again after their first description (Mammen, 1963;Watson, 2005). The other two wellestablished species, namely Cladocoryne floccosa Rotch, 1871 and Cladocoryne haddoni Kirkpatrick, 1890, can be easily distinguished by the number of aboral branched tentacles and the cnidome (Bouillon et al., 1987). The presence of unnamed Cladocoryne species has also been hypothesized by Schuchert (2003) while discussing Indonesian colonies described by Stechow and M€ uller (1923) and Vervoort (1941), yet no molecular analyses have been performed to assess the precise diversity of the genus. Moreover, the type material of C. floccosa, C. travancorensis, and C. littoralis cannot currently be located, and is presumably lost.
In contrast to Pteroclava Weill, 1931, the genus Cladocoryne is a generalist in terms of substrate choice, and has been reported on a variety of substrates, including rocks, algae, seagrasses, sponges, octocorals, bryozoans, polychaete tubes, and other hydroids (Bouillon et al., 1987;Gravili et al., 2015). Additional differences are found in the general morphology and reproductive structures, since Pteroclava has moniliform tentacles and produces free-swimming medusae, whereas Cladocoryne has oral capitate and aboral branched capitate tentacles and produces cryptomedusoids that do not detach from the parental polyp (Bouillon et al., 2006). However, both genera share the presence of euryteles grouped in rounded clusters in the polyp stage, a feature considered to be a synapomorphy of taxa ascribed to Cladocorynidae Allman, 1872 (Petersen, 1990). Other unique features of the family are (i) the presence of euryteles in the exumbrellar cnidocyst pouches, occurring in the freeswimming medusa of P. krempfi (Boero et al., 1995), and (ii) the branched aboral tentacles in all Cladocoryne species (Bouillon et al., 2006).
With this work we aimed at performing a thorough assessment of the family Cladocorynidae, including detailed morphological characterisations, phylogenetic analyses, single and multi-locus DNA-based species delimitations, statistics of morphological measurements, character evolution, and host-specificity assessment. Additionally, the phylogenetic position of the octocoral-associated species Zanclea timida was assessed based on new morphological information (i.e., the medusa stage) and genetic data.  Table S1). All species were collected in this depth range, with the only exception of Zanclea timida colonies at 40-50 m depth. When the presence of hydrozoan polyps was recorded in situ, small fragments of the substrate and associated hydroids were collected. Animals were anesthetized with menthol crystals, detached from their hosts using precision forceps, hypodermic syringe needles, and micropipettes and subsequently fixed in both 99% ethanol for molecular analyses and 10% formalin for morphological analyses. When possible, live animals were reared in constantly oxygenated bowls filled with seawater, at room temperature, under artificial light, and fed with Artemia nauplii until medusa release occurred.

Sampling and morphological assessment
Hydrozoans were identified to the species level following Boero et al., (1995), Bouillon et al., (1987), Schuchert (2006), and Puce et al., (2008). Cladocorynid specimens were observed and photographed under a Leica EZ4 D stereo microscope to assess their general morphology, and analyzed under a Zeiss Axioskop 40 compound microscope to study the fine-scale morphology of polyps, medusae, and cnidocysts. All measurements were taken using ImageJ 1.52p software. For the identification of host octocorals, small portions of tissue were immersed in a 10% sodium hypochlorite solution for 6 h to obtain sclerites. Octocorals were identified to genus level or, when possible, to the species level by examining the general morphology of the colony and the sclerites following Fabricius and Alderslade (2001) and Williams and Chen (2012).

DNA extraction and sequencing
Total genomic DNA was extracted from ethanol-fixed hydroids using the protocol described in Maggioni et al., (2020b). Additional DNA samples of Cladocoryne floccosa colonies from the Caribbean and Mediterranean Sea were obtained from Peter Schuchert (Mus eum d'Histoire Naturelle of Geneva, Switzerland, Table S1). Seven gene regions were amplified using the primers and protocols described in Maggioni et al., (2020c), namely portions of the mitochondrial large ribosomal RNA (16S rRNA), cytochrome oxidase subunit I (COX1), cytochrome oxidase subunit III (COX3), and nuclear small ribosomal RNA (18S rRNA), large ribosomal RNA (28S rRNA), internal transcribed spacer (ITS; including partial ITS1, 5.8S, and partial ITS2 regions), and histone H3 (H3) ( Table 1). PCR products were checked through 1.5% agarose electrophoretic runs, purified with Illustra ExoStar (GE Healthcare: Amersham, UK) and sequenced in both directions with ABI 3730xl DNA Analyzer (Applied Biosystems: Carlsbad, CA, USA). Geneious 6.1.6 was used to visually check, correct, and assemble the obtained chromatograms and to check for the presence of open reading frames in protein-coding genes. The consensus sequences obtained in this study were deposited in GenBank with the accession numbers listed in Table S1. Sequences of each DNA region were aligned using MAFFT 7.110 (Katoh and Standley, 2013) with the E-INS-i option, after adding outgroup sequences following Maggioni et al., (2017Maggioni et al., ( , 2018 (Tables S1). The ITS alignment was further run through Gblocks (Castresana, 2000;Talavera and Castresana, 2007) using the default 'less stringent' settings to remove ambiguously aligned regions. All alignments were also concatenated using Mesquite 3.2 (Maddison and Maddison, 2006).
Phylogenetic inference was performed using maximum parsimony (MP), maximum likelihood (ML), and Bayesian inference (BI) using both single-and multi-locus datasets and running the analyses on the CIPRES server (Miller et al., 2010). MP analyses were performed using PAUP 4.0b10 (Swofford, 2003) with heuristic searches stepwise addition and tree-bisection-reconnection branch swapping. Node support was assessed using 1000 bootstrap replicates with randomly added taxa. For ML analyses, the substitution model GTR+G was selected for each DNA region and partition, and analyses were run using RAxML 8.2.12 (Stamatakis, 2014) with 1000 non-parametric bootstrap replicates. BI analyses were performed using MrBayes 3.2.6 (Ronquist et al., 2012) using the substitution models and partitions listed in Table S2. Two independent runs for four Markov chains were conducted for 50 million generations, with trees sampled every 5000th generation. For all Bayesian analyses, parameter estimates and convergence were checked using Tracer 1.6 , where burn-in was set at 25%. MEGA X (Kumar et al., 2018) was used to calculate genetic distances within and among the obtained clades. Genetic distances were calculated as % uncorrected p-distances with 1000 non-parametric bootstrap replicates.
A species tree was obtained using *BEAST (Heled and Drummond, 2009) in BEAST 2.2.0 (Bouckaert et al., 2014). Specimens were assigned to different species according to the species delimitation analyses described below, and mitochondrial loci were considered as a single partition, whereas nuclear loci were kept separate. Yule process prior, together with a linear and constant-root population-size model, were used and each analysis was run for 10 8 generations, sampling every 10000th generation, and setting a 25% burn-in. The convergence of the analysis was checked using Tracer 1.6 , as done for the gene trees.
To better visualize the presence of host-related and geographyrelated genetic structure within Pteroclava clades, median-joining haplotype networks were built using the software PopART 1.7 (Leigh and Bryant, 2015). The most complete mitochondrial singlelocus dataset (16S rRNA) was used, with haplotypes coloured according to the genus of the host and the geographic provenance (Indian Ocean, Pacific Ocean, Red Sea).
Finally, a set of characters was mapped onto a reduced multilocus phylogenetic reconstruction of the Zancleida Russel, 1953 (including one specimen for each nominal species) obtained with BEAST 1.8.2 (Drummond et al., 2012) using the same parameters described for the *BEAST analysis. The mapped characters were: (i) the organisation of euryteles in the polyp stage (not organized, organized in clusters or a ring); (ii) the exumbrellar cnidocyst pouches in the medusa stage (absent, containing euryteles, containing stenoteles); (iii) the substrate choice of the polyp stage (generalist, symbiotic); (iv) the host of the polyp stage (none, octocorals, sponges, scleractinians, bryozoans). Stochastic mapping (Huelsenbeck et al., 2003) was used to map probable realizations of the evolution of the considered characters on the Zancleida tree. The analyses were run in the R environment (R Core Team, 2020), using the 'make.simmap' function in the package 'Phytools' (Revell, 2012). Reconstructions were performed under the 'equal rates' (ER) and the 'all rates different' (ARD) models, comparing the fit of the models with a likelihood ratio test using the function 'pchisq'. The ER model was selected for the characters 'organization of euryteles' and 'substrate choice', whereas the ARD model for the characters 'cnidocyst pouches' and 'host'. 10000 stochastic mapping replicates were conducted for each analysis and the results were summarized with pie charts representing the posterior probability of each internal node being in each state.

DNA-based species discovery and validation
Distance-and tree-based species delimitation approaches were used to investigate the presence of multiple species hypotheses in octocoral-associated Pteroclava. All analyses were run separately on the single-locus datasets, keeping only Pteroclava sequences and collapsing the alignments into haplotypes by removing identical sequences, using FaBox (Villesen, 2007). Sequences with different length but identical nucleotides were considered the same haplotype.
For tree-based methods, single-locus ML trees were obtained as described above and ultrametric Bayesian trees were obtained with BEAST 1.8.2 (Drummond et al., 2012), setting a coalescent tree prior and an uncorrelated lognormal relaxed clock. Three replicate analyses were run for 10 8 million generations, with trees sampled every 10000th generation, and were combined using LogCombiner 1.8.2 (Drummond et al., 2012) with a burn-in set to 25%. Maximum clade credibility trees were obtained using TreeAnnotator 1.8.2 (Drummond et al., 2012). The obtained trees were used as inputs for the Poisson Tree Process (PTP) and Generalized Mixed Yule Coalescence (GMYC) methods. Multiple-threshold PTP analyses (mtPTP: Kapli et al., 2017) were performed on the website 'mPTP Webservice' (https://mptp.h-its.org) using ML trees as input, whereas GMYC analyses were run in the R environment (R Core Team, 2020), using ultrametric trees as input. Specifically, single-threshold GMYC (stGMYC; Pons et al., 2006) and bGMYC (Reid and Carstens, 2012) analyses were run using the packages 'Splits' (Ezard et al., 2009), 'Ape' (Paradis et al., 2004), and 'bGMYC' (Reid and Carstens, 2012). bGMYC analyses were performed using a subset of 100 trees retrieved from the 10000 posterior trees obtained with each BI analysis. Single-threshold PTP and multiple-threshold GMYC methods were previously shown to be outperformed by the other implementations of PTP and GMYC used in this work (Fujisawa and Barraclough, 2013;Kapli et al., 2017) and were therefore not included in the analyses.
To validate the species hypotheses proposed by the single-locus analyses, the multi-locus dataset was analyzed under the multispecies coalescent model using Bayesian Phylogenetics and Phylogeography (BPP) 3.4 (Yang and Rannala, 2010). The algorithm A11 (joint species delimitation and species tree inference) was used and different values of the root age (s) and prior distributions of the ancestral population size (h) were tested, due to the lack of knowledge about these parameters in Pteroclava and to test their possible influence on the number of delimited species (Leach e and Fujita, 2010). We tested for: (i) large ancestral population size and deep divergence [h~IG (3, 0.04), s~IG (3, 0.02)]; (ii) small ancestral population size and shallow divergence [h~IG (3, 0.004), s~IG (3, 0.002)]; (iii) small ancestral population size and deep divergence [h~IG (3, 0.004), sĨ G (3, 0.02)]; (iv) large ancestral population size and shallow divergence [h~IG (3, 0.04), s~IG (3, 0.002)]. Each analysis was run using both the reversible-jump MCMC algorithm 0 (ɛ = 1) and algorithm 1 (a = 1.5, m = 1.5) to assess the impact of the fine-tuning parameters of reversible-jump algorithms on the speciation probabilities, and was repeated two times, to confirm the stability of results across runs. For each run, 50000 samples were collected with a sampling frequency of 10 iterations after a burn-in of the first 5000 cycles.
Finally, DNA diagnostic characters were searched for each P. krempfi clade, using the package QUIDDICH (K€ uhn and Haase, 2020) in the R environment (R Core Team, 2020). Specifically, character attributes (i.e., single nucleotides) present in all members of a defined clade but absent in members of other clades were searched, and the analyses were run for each molecular marker, considering different species hypotheses.
Alignments and raw phylogenetic trees are provided as supplementary files.

Pteroclava cnidocyst size analyses
To explore differences in the size of cnidocysts among P. krempfi clades, we selected 15 colonies, five for each of the clades Ia, IIa, and III, whereas measurements were not performed for clades Ib and IIb (Fig. 2), since insufficient material was available. The analyzed colonies were associated with nine different hosts and collected from five localities. For each colony, ten polyps were selected and, for each polyp, up to 15 capsules per type were measured. The descriptive statistical parameters of all cnidocyst types [length, width (mean AE sd), range (maximum-minimum)] were calculated for each clade in the R environment (R Core Team, 2020) using the 'Rcmdr' package (Fox and Bouchet-Valat, 2020). Kernel density plots (Sheather, 2004) were produced for the length of each cnidocyst type and clade. Stratified box plots for each cnidocyst type, clade, locality and host were also obtained. All graphics were made using the package 'ggplot2' (Wickham, 2016).
To study the variation of cnidocyst sizes in the P. krempfi clades, linear mixed models (LMM) or generalized linear mixed models (GLMM) were fitted for each of the three cnidocyst types (eurytele, large stenotele, small stenotele). Only the length data were used, since width data are in general less variable (Garese et al., 2016). The selection of LMM or GLMM was performed assessing the normality or not of the residuals of the models (following Garese et al., 2016). The models were performed using the R package 'lme4' (Bates et al., 2015). The general form of the models was: where y is the response variable, b i is the fixed effect parameter, µ i is the random effect parameter, x is the matrix of fixed parameters, z is the matrix of random effects and e i is the matrix of errors. The global model form included 'cnidocyst size' as the dependent variable and 'clade' as fixed effect. Since datasets were assembled by measuring polyps from different colonies, localities, and/or hosts, these variables were at first considered as random effects and then their significance evaluated. By consequence, the global model form was: Cnidocyst sizes~Clade + (1|Polyp) + (1|Colony) + (1| Locality) + (1|Host) + e.
To find the best model, the significance of each random effect was tested by deleting them one by one, using 'ranova' function of the R package 'lmerTest' (Kuznetsova et al., 2017), and to rank LMMs or GLMMs the Akaike Information Criterion (AIC) was used. Once the best model was determined, the normality of residuals of the LMM were tested by Shapiro-Wilks test and, graphically, in QQplots, and the homoscedasticity of variance was checked in fitted vs residual graphs. When the assumption of normality was not met, a GLMM was fitted by maximum likelihood (Laplace Approximation) with Gamma distribution for errors and inverse link function. Then, the best models were compared through ANOVA with a null model (Cnidocyst sizes~Clade) to test the significance of the LMM or GLMM. Finally, the model estimates for each clade, the confidence intervals, and the standard deviation for random effects, were extracted and analyzed. The R scripts used for the LMM and GLMM are shown in File S1.

Molecular phylogenetics and species delimitation
The length of single-locus alignments was 625 bp for the 16S rRNA, 677 bp for COX1, 668 bp for COX3, 1712 bp for 18S rRNA, 1684 bp for 28S rRNA, 662 bp for ITS (927 bp before the Gblocks treatment), and 352 bp for H3 (Table 1), resulting in a total of 6380 bp for the concatenated alignment. Overall, the molecular markers were amplified and sequenced with high success (Table 3) and approximately 80% of the terminals in the concatenated dataset were represented by all markers. When this was not the case, one to three markers were missing (Table S1). The single-and multi-locus phylogeny reconstructions were broadly concordant (Fig. 2a, b; Fig. S1; Table S3). Although the relationships among clades varied across the single-locus analyses, each of the analyzed colonies belonged to the same molecular clade in almost all trees and phylogeny reconstruction criteria (Table S3). Multi-locus analyses resulted in the same ingroup topology with high nodal support (Fig. 2a, b), with only minor variation in the relationships among outgroups ( Fig. 2a; Fig. S1). Here we considered the node support as maximal when bootstrap values (BS) for MP and ML analyses were = 100 and Bayesian posterior probabilities (BPP) = 1, high when BS ≥ 75 and BPP ≥ 0.9, and low when BS < 75 and BPP < 0.9. The family Cladocorynidae was monophyletic in all multilocus analyses, with maximal support in ML and BI analyses, and lower support in the MP analysis (BS = 71), and its phylogenetic position as sister group of the families Asyncorynidae Kramp, 1949, Milleporidae Fleming, 1828, Solanderiidae Marshall, 1892, and Zancleidae Russel, 1953 was confirmed (Fig. 2a,  Fig. S1). Both cladocorynid genera, Cladocoryne and Pteroclava, and all the analyzed species, were consistently recovered as monophyletic with high or maximal support values (Fig. 2b). All the specimens belonging to Pseudozanclea timida, and collected in the Maldives, formed a monophyletic clade that was sister to Pteroclava in all multi-locus analyses with maximal nodal support. The phylogenetic placement of this species within the Cladocorynidae was therefore demonstrated, supporting the transfer from the genus Zanclea and family Zancleidae to the new genus Pseudozanclea. In all analyses, Pteroclava krempfi was characterized by a strong intra-specific genetic variation, spanning in most cases five well-supported clades ( Fig. 2b; Fig. S1; Table S3): Pteroclava krempfi Ia, Ib, IIa, IIb, and III. The clades P. krempfi I, II, and III were recovered in all analyses, whereas sub-clades Ia, Ib, IIa, and IIb were recovered in most analyses, with the exception of single-locus ITS, 28S, and H3 analyses (Table S3).
The genetic distances among the nominal species were high for all DNA regions, with nuclear 18S and 28S rRNA showing the lowest values, whereas intraspecific distances were low for all species but P. krempfi and C. floccosa (Table 4; Table S4). These two species showed high intra-specific distances especially in mitochondrial regions, with values of 4.7% for P. krempfi and 3.7% for C. floccosa for the 16S rRNA dataset (Table S4). Regarding P. krempfi, the inter-clade genetic distances were also generally high, Bayesian posterior probabilities of BI analyses, respectively, for the multi-locus analyses, with asterisks indicating maximal nodal support for all analyses; nodal supports obtained with single-locus analyses are shown in panels, coloured as coded in the legend. In (c) numbers only represent Bayesian posterior probabilities. In (b) clades are highlighted with different colours and alphanumeric codes; sampling localities, as coded in the legend, and hosts are also reported for each terminal taxon. and a reasonable decrease in the intra-clade distances was observed by grouping sequences in the five clades obtained with phylogenetic analyses, with no overlap with intra-clade distances (Table 4; Table S4).
Single-locus DNA-based species delimitation analyses of P. krempfi showed variable results, according to the method used and DNA region analyzed (Fig. 3a), and resulted in five possible scenarios, spanning one to five species hypotheses (Fig. 3b). According to phylogenetic and genetic distance analyses, the delimitation based on mitochondrial regions revealed the highest number of species hypotheses, whereas those based on nuclear regions mostly retrieved three species hypotheses (Fig. 3a). The BPP multi-locus delimitation analyses consistently inferred a five-species model with high posterior support (always > 0.99), regardless of the algorithm and sets of priors used. Consistently, the posterior probabilities for each of the five delimited species hypotheses were always high (> 0.99) (Table S5).
Diagnostic characters were searched in P. krempfi clades considering the three-species model and the fivespecies model (Table S6). In the three-species model, diagnostic characters were detected for all clades and molecular markers, with highest numbers in the mitochondrial COX1 and COX3 genes and the ITS. On the other hand, in the five-species model diagnostic characters were found for all clades only in the 16S rRNA and COX3 genes. Overall, clade III consistently showed the highest number of diagnostic characters in all DNA regions.
A species tree of the Cladocorynidae family was produced considering P. krempfi as composed of five  The column 'Conc' refers to the % of terminals with all markers sequenced in the concatenated dataset species (Fig. 2c). The topology of the species tree was identical to that of the concatenated multi-locus trees, although an overall lower support was obtained (Fig. 2b, c). Specifically, the relationships among Pteroclava clades did not vary, Pseudozanclea timida was confirmed as the sister group of Pteroclava, and Cladocoryne as the sister group of Pteroclava + Pseudozanclea.
According to the ancestral state reconstructions, the organisation of euryteles in clusters or rings in the polyp stage (Fig. 4a) and the presence of exumbrellar cnidocyst pouches with euryteles in the medusa (Fig. 4b) resulted to likely be synapomorphies of the Cladocorynidae, being appeared in the most recent common ancestor (MRCA) of the extant family members. The establishment of symbiotic relationships between the polyp stage and other benthic organisms used as substrates emerged multiple times in the superfamily Zancleida (Fig. 4c), but the association with octocorals probably emerged in the MRCA of the Pteroclava + Pseudozanclea clade (Fig. 4d).
Distribution, host specificity, and geographic structure of Pteroclava clades Both Pteroclava krempfi clades I and II were found across the Indo-Pacific and Red Sea, whereas clade III was exclusively observed in Caribbean localities and associated with Antillogorgia. Clade I was associated with members of the octocoral family Alcyoniidae (Sinularia, Sarcophyton, Lobophytum, Rhytisma), with the only exception of a colony associated with Paralemnalia (family Nephtheidae) ( Table 2). Clade Ib included four colonies collected in the Red Sea and Australia, associated with Sarcophyton and Sinularia octocorals, whereas the remaining specimens belonged to clade Ia. Within clade I, no clear host-related genetic structure was observed (Fig 4a), with two haplotypes associated with multiple hosts, whereas populations from the Indian Ocean and Red Sea appeared to be moderately differentiated, with the exception of a Red Sea haplotype more related to Indian Ocean haplotypes (Fig. 5b). Clade II was composed of clade IIa, associated with the genera Astrogorgia and Paraplexaura (Plexauridae), Melithaea (Melithaeidae), Ctenocella (Ellisellidae) and Dendronephthya (Nephtheidae), and clade IIb, represented by a single colony associated with Acanthogorgia (Acanthogorgiidae) ( Table 2). In clade II no host-related genetic structure was also clearly detectable, and the most common haplotype was associated with two different octocoral genera (Fig. 5c). Clade II showed a lower number of haplotypes compared to clade I and in this case different localities did not share the same haplotypes, even if a clear geographic pattern was not observed (Fig. 5d).

Morphometry of the Pteroclava polyp cnidome
The cnidome of P. krempfi polyps was composed of stenoteles of two size classes (small and large) located in the tentacles and macrobasic apotrichous euryteles below the hypostome. A total of 5816 capsules were measured (Table S7), and the raw data reflected similarities in cnidocyst size among the three investigated clades (Ia, IIa, and III), as shown in Table S8. The eurytele length of clade III was slightly larger than Ia and IIa, even if the ranges of the three clades were partially overlapping. Furthermore, eurytele width was slightly larger in clade III. All other measurements were not noticeably different among the investigated clades.
The distributions of cnidocyst length for each cnidocyst type and clade are reported in the density plots in Fig. 5a-c, showing the general overlapping of the length distributions and the higher values for eurytele length in clade III (Fig. 6a). Fig. 5d-f show the stratified box plots of the length measurements for each cnidocyst type, clade, host, and locality investigated. The larger size of euryteles of clade III can be observed in Fig. 6d, whereas clade Ia from the Red Sea showed the smallest size. The size of both large and small stenoteles appeared to be very similar  among the three clades from different hosts and localities (Fig. 5e, f), with the largest values observed in specimens belonging to clade IIa associated with Melithaea in the Red Sea. The best model produced by model selection was: Cnidocyst sizes~Clade + (1|Polyp) + (1|Colony) + e. For all cnidocyst types, both 'Polyp' and 'Colony' were significant, whereas 'Locality' and 'Host' were not (Table S9). For euryteles and large stenoteles, a LMM was fitted according to the normality of the residuals (P = 0.4805 and P = 0.2207, respectively; a = 0.05) (Fig. S2) from the best models. For small stenoteles, the residuals of the best LMM did not fit the normal distribution (P = 0.001; a = 0.05), therefore a GLMM was fitted. The best models for eurytele and stenotele large (LMMs), and stenotele small (GLMM) were significant versus the null models (Table 5). Estimates and confidence intervals (CIs) of the models for each clade are shown in Table 6, and the standard deviation of the random effects 'Polyp' and 'Colony' for each cnidocyst type are shown in Table S10.
Euryteles of the clade III showed an estimated value for cnidocyst length larger than Ia and IIa. Moreover, the CIs for clade III did not overlap with the CIs of clade Ia, and only partially overlapped with the CIs of clade IIb. On the other hand, estimates for clade Ia and Ib were very similar, with almost totally overlapping CIs (Table 6). Although both 'Polyp' and 'Colony' as random effects were significant, their standard deviations were small, around 1 µm and 2 µm, respectively (Table S10), suggesting no considerable variation within polyps and colonies. The model estimates for both small and large stenoteles were almost identical for the three clades, with completely overlapping CIs and very small standard deviations of 'Polyp' and 'Colony' (Table 6; Table S10). Thus, no differences were detected in the length of stenoteles among clades.
Description: Colonies living as partial endosymbionts of a variety of octocorals (Table 2, Fig. S3). Unbranched stems (up to 1.5 mm long) arising from creeping hydrorhiza embedded in host tissues (Fig. 9a). Perisarc stopping at the base of polyps (Fig. 9b). Polyps claviform up to 2 mm high, with 4-5 oral (up to 0.5 mm long) and 6-24 aboral (up to 1 mm long) moniliform tentacles scattered on the hydranth (Fig. 8a, b) with stenoteles of two size classes. Up to four rounded patches of macrobasic apotrichous euryteles below the oral tentacles, sometimes absent. Up to five medusa buds at the base or half of the polyp body, polyps often showing reproductive exhaustion (Fig. 9c). Newly liberated medusa (Fig. 9d) with a bell-shaped umbrella (up to 800 lm high and 850 lm wide), manubrium extending for 1/3 or 1/2 of the umbrella, and microbasic euryteles scattered on the exumbrella. Four perradial canals ending in two tentaculate bulbs and two smaller atentaculate bulbs. Atentaculate bulbs with an exumbrellar cnidocyst pouch containing up to four euryteles (Fig. 9e), the latter often found also in tentaculate bulbs. Tentacular bulbs bearing one tentacle each (Fig. 9f), about 1 mm long, with up to 40 ciliated rounded cnidophores approximately 25 lm wide, armed with 4-5 bean-shaped macrobasic apotrichous euryteles. Living polyps transparent, with a whitish hypostome and sometimes a reddish gastroderm (Fig. 9a), medusae transparent, with whitish or reddish manubrium.
Remarks: The establishment of the new genus is supported by both morphological and genetic data. Indeed, the polyp stage of Pseudozanclea differs from Pteroclava and Cladocoryne by showing a Zanclea-like morphology (i.e., hydranths with oral and aboral capitate tentacles) and the presence of peculiar tissue pockets containing medusa buds. Similarities with Pteroclava are observed in the newly released medusa, with exumbrellar pouches containing euryteles, and in the association with octocorals. These latter similarities are also reflected in the proposed phylogeny of the Cladocorynidae, with a sister group relationship between Pteroclava and Pseudozanclea.
Description: Maldivian colonies associated with the octocoral Bebryce cf. grandicalyx and with a sponge in the family Suberitidae (Fig. 10a). Indonesian colonies associated with the octocoral Paratelesto sp. Perisarcfree hydrorhiza with clusters of macrobasic apotrichous euryteles (Fig. 10b), extending on the octocoral surface and completely covered by the sponge, when present. Polyps piriform, up to 0.9 mm long, highly contractile (Fig. 9b, c), and partially covered by the sponge (Fig. 10d). Five or six oral capitate tentacles up to 200 lm long (capitula~60 lm wide) and 9-12 aboral tentacles up to 180 lm long and with smaller capitula (capitula~45 lm wide). Capitula containing stenoteles of two size classes. Polyps with a proximal ring of macrobasic apotrichous euryteles (Fig. 10e) corresponding to the point where the sponge stops surrounding the polyps, approximately at one-third of the polyp body. Up to four medusa buds borne basally, kept inside a tissue pocket protected externally by the ring of euryteles (Fig. 10e) and covered by the sponge. Medusa buds exposed only at maturation, before release. Newly liberated medusa (Fig. 9f, g) with a bell-shaped umbrella (~750 lm high and~800 lm wide), manubrium extending for one-third or half of the umbrella (Fig. 10g), and microbasic euryteles scattered on the exumbrella. Four perradial canals ending in two triangular tentaculate bulbs and two smaller atentaculate bulbs (Fig. 9g, h). Atentaculate bulbs with an exumbrellar cnidocyst pouch containing up to five euryteles, and in some cases small-sized stenoteles (Fig. 10i). Euryteles found also in tentaculate bulbs (Fig. 10h). Tentacular bubs bearing one tentacle each, about 700 lm long, with 50-100 rounded ciliated cnidophores approximately 30 lm wide, armed with 3-4 bean-shaped macrobasic apotrichous euryteles. Seven days old medusa with similar size and longer tentacles (up to 2 mm) equipped with more cnidophores. Living polyps transparent (Fig. 9c, d), sometimes with a reddish gastroderm and living medusae transparent, with reddish bulbs and a whitish manubrium (Fig. 10f).

An updated phylogeny of the Cladocorynidae
In this work, a phylogenetic hypothesis of the family Cladocorynidae was proposed and, according to multiple lines of evidence, the new genus Pseudozanclea was established to accommodate Zanclea timida. Various phylogenetic analyses consistently placed Pseudozanclea timida as the sister group of Pteroclava and as phylogenetically distant from other Zanclea species. The description, for the first time, of the newly released medusa of P. timida demonstrated that there are strong similarities between the Pseudozanclea and Pteroclava medusae, both in the general morphology and in the unique presence of exumbrellar cnidocyst pouches containing euryteles. The polyp stage also shares a similarity with those of other cladocorynid taxa, which is the presence of macrobasic apotrichous euryteles organized in a specific pattern. Finally, similarly to Pteroclava, Pseudozanclea is an obligate symbiont of octocorals. The family Cladocorynidae was initially erected to accommodate the genus Cladocoryne (Allman, 1872), whereas only later Petersen (1990) and Boero et al., (1995) proposed to include the genus Pteroclava. This conclusion was justified by the synapomorphy represented by the presence of patches of large macrobasic euryteles in the polyp stages of both genera. However, similar patches are found at least in another noncladocorynid species, Eudendrium glomeratum Picard, 1952, in which patches can form an irregular band (Schuchert, 2008). The presence of eurytele patches in Cladocorynidae and E. glomeratum clearly represents morphological convergence, with the trait emerging twice independently in two phylogenetically distant suborders (Maronna et al., 2016). As demonstrated in this work, Pseudozanclea timida shows aggregations of macrobasic euryteles in the polyp stage, even if organized in a ring rather than patches, and this organization may derive from the coalescence of the patches observed in Pteroclava and Cladocoryne. Interestingly, specimens from Indonesia and associated with Paratelesto possess macrobasic apotrichous mastigophores instead of euryteles, even if the heteronemes of the two populations show similarity in shape and overlap in size. Maldivian heteronemes have an enlargement of the distal part of the shaft, typical of euryteles, and similar to that occurring in Pteroclava and Cladocoryne, i.e., with a zig-zag pattern and covered with spines (see for instance Maggioni et al., 2016: 488, Fig. 7 andBouillon et al., 1987: 283: Fig. 7). On the other hand, the shaft figured by Puce et al. (2008b: 1652) appears to be the same diameter for all its length, a typical feature of mastigophores. This could be the case, but it is noteworthy that the figured distal portion of the shaft seems damaged, for instance with spines not covering the end of the shaft, and this may have caused an incorrect classification of the cnidocyst occurring in Indonesian specimens.
Whatever the cnidocyst type found in the ring, the function may relate to the polyp defence, since P. timida hydroids are able to retract within their basal cup, exposing only the cnidocyst ring and the tentacular capitula (Puce et al., 2008b). Moreover, the immature medusa buds are kept inside a tissue pocket below the cnidocyst ring, and the latter may provide protection to developing medusae. Pseudozanclea timida has stolonal cnidocyst clusters, which have been hypothesized to protect the perisarc-derived stolons from predators (Puce et al., 2008b). However, the stolons of Maldivian specimens were always covered by a sponge, contrary to those described from Indonesia. Therefore, the defensive function of stolonal cnidocyst clusters does not seem to apply for the Maldivian specimens of P. timida. The associated sponge also covers the proximal portion of P. timida hydranths, including the medusa buds, and this may result in a further protection of the immature medusae before release (Maggioni et al., 2020a). The sponge was never observed to cover P. timida polyps above the cnidocyst ring, and the aggregation of euryteles may prevent an excessive overgrowth of the hydroids by the sponge.
Another synapomorphy of the family is the presence of exumbrellar pouches with eurytele capsules, a peculiarity never found in other hydrozoan species. The medusae of both Pseudozanclea timida and Pteroclava krempfi clearly possess this feature on atentaculate bulbs. Cladocoryne species also show exumbrellar eurytele aggregations on their cryptomedusoids. The presence of these exumbrellar cnidocyst clusters was confirmed herein for both C. floccosa and C. haddoni, and are apparently present also in C. minuta (Watson, 2005), even if described as scattered rather than grouped. Regarding the other two Cladocoryne species, the reproductive structures of C. travancorensis are unknown, whereas the cnidome of the C. littoralis cryptomedusoids was not reported in the original description (Mammen, 1963). It is therefore likely that the exumbrellar cnidocyst pouches containing euryteles are found in all cladocorynid species, but future reanalysis of the poorly described species will be necessary to address this issue.
The clade Pseudozanclea + Pteroclava is composed of obligate octocoral symbionts, and this association is likely to have arisen in the MRCA of the two genera. The species Pteroclava crassa was described growing on the hydroid Macrorhynchia philippina Kirchenpauer, 1872 (Pictet, 1893), but the information on this species is very scant and the only difference with P. krempfi appears to be the host. Therefore, it remains unclear whether the two Pteroclava species are conspecific or not. Within the superfamily Zancleida, another species was described growing on octocorals, Zanclea cubensis Varela, 2012, described from the Caribbean Sea (Varela, 2012). Zanclea cubensis is very similar to other zancleid species, having a claviform hydranth with oral and aboral capitate tentacles, and stenoteles and euryteles as cnidocysts. The medusa was not described, and the species was named mostly based on the combination of the type of euryteles (macrobasic holotrichous) and the association with octocorals. Considering the lack of genetic data and the partial description of the species, along with the polyphyletic nature of the family Zancleidae and the genus Zanclea , the phylogenetic position of Z. cubensis remains unknown. Therefore, it is currently not possible to establish how many times the association with octocorals evolved in the Zancleida.

Taxonomic uncertainties in the Cladocorynidae
Previous studies suggested the presence of three cryptic species within Pteroclava krempfi, based on the analysis of four DNA regions of Maldivian and Caribbean samples Montano et al., 2017). In the present work, several colonies associated with multiple hosts and from different localities were analyzed to better comprehend the genetic structure of the P. krempfi species complex. Despite the conserved general morphology across samples, our phylogenetic, genetic distance, and single-and multi-locus species delimitation analyses confirmed the presence of highly divergent genetic lineages. The number of possible species hypotheses varied according to the method used and DNA region analyzed, with mitochondrial and multi-locus analyses resulting in the highest numbers of species hypotheses. The higher number of species found based on mitochondrial regions is concordant with the higher rates of nucleotide substitution found in mitochondrial DNA of hydrozoans compared to other phylogenetically informative nuclear regions (e.g., Maggioni et al., 2020b, c). It therefore remains unclear how many possible species exist under the name P. krempfi, also because clades Ib and IIb are herein represented by a limited number of samples. Unfortunately, it is currently not possible to formally describe the clades due to the lack of P. krempfi material suitable for genetic analyses from the type locality (Vietnam) and associated with the type's host (Cladiella krempfi). Additionally, Cladocoryne floccosa may be another species complex, given the high intraspecific genetic distances found in this study for all mitochondrial and some nuclear regions. Cladocoryne floccosa has a wide distribution, from tropical to temperate areas (Schuchert, 2006;Gravili et al., 2015), and only a thorough sampling across the whole geographic range and an accurate morpho-molecular assessment would clarify the possible presence of cryptic lineages or species.
The analysis of cnidocyst size data of P. krempfi revealed partial differences among the genetic clades, since a clear size variation was found only for euryteles, showing a larger size in organisms from the Caribbean (clade III). However, no differences were detected in eurytele size between clade Ia and IIa and stenotele size among all clades. The statistical treatment of cnidocyst size data was previously used to search for fine-scale variations in other hydrozoan species, revealing in some cases no differences (Wollschlager et al., 2013) or significant differences among genetic clades Manca et al., 2019;Maggioni et al., 2020c). Similar studies were performed also on sea anemones, in which only some cnidocysts, coming from particular structures, varied among morphs (Gonz alez-Muñoz et al., 2017), or did not vary at all (Gonz alez-Muñoz et al., 2018), showing a high intraspecific variability (Garese et al., 2016). The results obtained in this work are therefore in line with previous works showing limited, but in some cases still useful, power of discrimination of cnidocyst size among closely related species or populations. Nevertheless, to present, there is no evidence that can explain the intraspecific variation on cnidocysts.
The obtained P. krempfi clades showed different host preferences. Clades Ia and IIb showed overlapping host octocoral genera, and all clade I colonies but one were found living on genera belonging to the family Alcyoniidae. The only exception was a colony found on Paralemnalia (family Nephtheidae). No host overlapping was found between clade I and II at the genus level, with a single exception at the family level being a clade IIa colony associated with Dendronephthya, another genus of the family Nephtheidae. All other clade II colonies were associated with several octocoral genera and families, and no host overlapping between clades IIa and IIb was observed at both the genus and family level, even if clade IIb is currently composed of a single colony. Finally, clade III was consistently found associated with Antillogorgia ( Montano et al., 2017;Miglietta et al., 2018). Another P. krempfi record from the Caribbean Sea reports a colony associated with Plexaurella grisea Kunze, 1916 (Varela andCabrales Caballero, 2010), suggesting a wider host range for clade III. Overall, a certain degree of host specificity can be detected, mostly at the family level, even if with some exceptions. Host specificity has recently been demonstrated for several other coral-reef invertebrates associated with a variety of benthic organisms, such as parasitic gastropods (Gittenberger and Gittenberger, 2011;Gittenberger and Hoeksema, 2013;Potkamp et al., 2017;Fritts-Penniman et al., 2020), coral-dwelling barnacles (Malay and Michonneau, 2014;Tsang et al., 2014), commensal shrimps (Hork a et al., 2016), copepods (Korzhavina et al., 2019), benthic ctenophores (Alamaru et al., 2017), coral gall-crabs (Garc ıa-Hern andez et al., 2020), acoel flatworms (Kunihiro et al., 2019), and hydrozoans (Maggioni et al., 2020b, c). Hostspecificity of associated species can be very weak as observed in various coral-dwelling copepods living on mushroom corals (Ivanenko et al., 2018) and in one serpulid worm living on a wide range of Caribbean scleractinians . Hence, the causes (and consequences) of these, more or less specific, associations largely remain unclear and need further research. Additionally, colonies belonging to the same clade (clade I) displayed variable host preference, with some Alcyoniidae genera showing higher prevalence of the symbiosis, and these prevalences also varied according to the studied locality and reef type Seveso et al., 2020). Therefore, despite a single clade being associated with multiple octocoral genera, it is likely that, at least in some localities, some of the associations are rarer than others, or even occasional, and that the association is host-reliant .
In conclusion, this work sheds new light on the evolution and diversity of the poorly known family Cladocorynidae, providing a new well-supported phylogenetic reconstruction of the family, establishing a new genus, and demonstrating the presence of several DNA-based species hypotheses. Moreover, the new data on the host specificity and geographic distribution of P. krempfi provide valuable insights for future research on how the current diversity patterns have evolved and on the study of the nature of the hydrozoan-octocoral symbioses.

Supporting Information
Additional supporting information may be found online in the Supporting Information section at the end of the article. Fig. S1. Phylogentic reconstructions based on the multi-locus dataset according to (a) Maximum parsimony, (b) Maximum likelihood, (c) Bayesian inference. Fig. S2. Graphical evaluation of normality of residual to the fit of a linear model. QQ-plots (left) and dispersion diagram plots (right) for each cnidocyst type. Fig. S3. Overview of the association between Pteroclava krempfi and alcyonancean octocorals. Table S1. Table S1. Information on the specimens included in the analyses. Newly obtained sequences are in bold. Table S2. Partitions and substitution models based on AIC for the multi-locus dataset. Table S3. Support values for the clades recovered by single-and multi-locus analyses, based on maximum parsimony (MP), maximum likelihood (ML), Bayesian inference (BI). Highly supported nodes are in green (MP and ML bootstrap values ≥ 75, BI posterior probabilities ≥ 0.9). Table S4. Genetic distances shown as % uncorrected p-distance (AE standard deviation) among and within cladocorynid species/clades. Table S5. Outputs of the A11 BPP analyses based on different sets of priors and algorithms. Posterior probabilities (PP) are shown when > 0. Table S6. Molecular diagnostic characters for each Pteroclava krempfi clade in the 3-and 5-species hypotheses for (a) 16S rRNA, (b) COX1, (c) COX3, (d) 18S rRNA, (e) 28S rRNA, (f) ITS, and (g) H3 markers. n.s. not sequenced. Table S7. Measurements of cnidocyst length and width in Pteroclava krempfi. Table S8. Descriptive statistical parameters of each Pteroclava krempfi cnidocyst type for clades Ia, IIa, III reported in µm as range (mean AE SD). Table S9. Evaluation of the significance of the random effects by deleting one by one from the complete model. Table S10. Pteroclava krempfi cnidocyst standard deviation (SD) of random effects from the best models.
File S1. R scripts used for the LMM and GLMM.