An integrative approach assesses the intraspecific variations of Procamallanus (Spirocamallanus) inopinatus, a common parasite in Neotropical freshwater fishes, and the phylogenetic patterns of Camallanidae

Abstract Integrative taxonomy was used to evaluate two component populations of Procamallanus (Spirocamallanus) inopinatus in Brazil and the phylogeny Camallanidae. Parasite populations were collected in the characiform Anostomoides passionis from River Xingu (Amazon basin) and Megaleporinus elongatus from River Miranda (Paraguay basin). Morphology was analysed using light and scanning electron microscopy (SEM). Genetic characterization was based on partial sequences of the 18S and 28S rDNA, and COI mtDNA. Phylogenies were based on 18S and COI due to data availability. Generalized Mixed Yule Coalescent (GMYC), Poisson Tree Process (PTP) and *BEAST were used for species delimitation and validation. SEM revealed for the first time the presence of minute denticles and pore-like structures surrounding the oral opening, phasmids in females and confirmed other important morphological aspects. Statistical comparison between the two-component populations indicated morphometric variations, especially among males. The different component population of P. (S.) inopinatus showed variable morphometry, but uniform morphology and were validated as conspecific by the GMYC, PTP and *BEAST. Some camallanid sequences in GenBank have incorrect taxonomic labelling. Host, environment and geographic aspects seem to be related to some lineages within Camallanidae; however, their real phylogenetic meanings are still unclear.


Morphological procedures
Nematodes were cleared in glycerine and observed in a Leica DM5500 B microscope, with DIC and the software LAS Leica TM , for morphometric analysis. Measurements are given in micrometres, unless otherwise stated. Specimens used for SEM (1 male and 1 female of each infrapopulation) were dehydrated through a graded ethanol series, dried by evaporation with hexamethyldisilazane, coated with gold and observed in a JEOL JSM 6460-LV, at an accelerating voltage of 15 kV. The systematic classification of camallanid nematodes follows the proposal of Moravec and Thatcher (1997). Voucher specimens were deposited in the Coleção Helmintológica do Instituto Oswaldo Cruz (accession nos. CHIOC 38939, 38940).

Statistical analysis of morphological data
In order to evaluate significant differences among morphometric features, we performed a series of statistical analyses. First, one sample Kolmogorov-Smirnov (K-S) test was performed to verify the normal distribution of the quantitative variables (Zar, 2010).
The quantitative variables were the standard length (SL) of fishes, and total body length, maximum width, buccal capsule length and width, distance from anterior end to nerve ring, deirids and excretory pore, total length of oesophagus, tail length, distance from cephalic end to vulva in females and spicules in males of parasites. All these measurements were also used for calculating their ratio to total body length; the coefficient of variation (CV) for each of the previous mentioned features was also calculated. These morphometric features were chosen because they are used in the specific diagnosis, and ratios were used for verifying if the organ's anatomy (or their locations on the cuticle surface) undergoes modifications according to the parasite growth. Ratios and CV are expressed in percentage (%) and mean values are followed by ± 1 standard deviation, unless otherwise stated.
According to the results of K-S test, we performed non-parametric inferential methods and generalised linear models (GLM). As a preliminary analysis, in order to verify possible influences of host-environment attributes and infection burden on the parasite body length, host SL and parasite intensity of infection (see Bush et al. 1997) were tested as predictive of parasite body length using linear regression (Zar, 2010).
First, host SL and parasite intensity of infection were tested separately against parasite body length using simple linear regression, in order to evaluate which variable was most predictive through the r 2 (coefficient of determination). Then a multiple linear regression was performed, inserting the predictive variables into the expected predictive model, in which the most predictive was inserted first followed by the less predictive (Zar, 2010). Following the previous methodology, but using logistic regression instead of linear, we tested the effect of host species and parasite sexual maturation (i.e. male, non-gravid and gravid females) on parasite body length, estimating the odds ratio (OR, 0 < OR < 1 indicates antagonistic effect; OR = 1 lack of effect; OR > 1 synergetic effect) and the confidence interval (CI) of each model (Dohoo et al. 2003). We also performed a Mann-Whitney (U) to evaluate and prove the differences in SL among the two fish species (Zar, 2010).
Following the preliminary analysis, parasites were divided into three categories according to sexual maturation: males, non-gravid females (without larvae in uterus) and gravid females (uterus filled with larvae). All morphometric features previously stated for parasites were tested against their total body length using the Spearman's https://www.cambridge.org/core/terms. https://doi.org/10.1017/S0031182020001687 Downloaded from https://www.cambridge.org/core. IP address: 201.231.191.128, on 16 Sep 2020 at 14:42:30, subject to the Cambridge Core terms of use, available at correlation (r s ) (Zar, 2010), considering all sexual maturation categories. Differences in morphometric features of parasites among the two component populations, i.e. from A. passionis and M. elongatus, were evaluated using the U test, within each sexual maturation category and among non-gravid and gravid females (Zar, 2010). Additional terminology related to parasite ecology follows the proposal of Bush et al. (1997).  (Chang et al. 2014). Saturation of nucleotide substitution was tested using Xia's method implemented in DAMBE (Xia et al. 2003;Xia, 2018).

Genetic procedures
Datasets were subjected to maximum likelihood (ML) and Bayesian inference (BI) analyses, using PHYML and MrBayes, respectively (Huelsenbeck and Ronquist, 2001;Guindon and Gascuel, 2003 generations and discarding the initial ¼ of sampled trees (1 × 10 6 ) as bur-in. In order to check chain convergence, analyses were run in duplicates and inspected using Tracer (Rambaut et al. 2018).
The genetic markers profile of phylogenetic informativeness was evaluated by PhyDesign (Townsend, 2007) based on the concatenated dataset (18S + COI). For this analysis, an ultrametric tree was generated using BEAST 2.5 (Bouckaert et al. 2019), the best-fit substitution model was chosen according to bModelTest (Bouckaert and Drummond, 2017), models of each partition were unlinked, clock model was selected based on marginal likelihood estimated (MLE) from stepping stone and path sampling (Baele et al. 2012) being 50 path steps, 5 × 10 5 itinerations and sampling every 500 generations. MCMC chains were run for 5 × 10 7 generations, sampling every 5 × 10 3 generations and ¼ of initial runs discarded as burn-in. All analyses were run in duplicates and inspected with Tracer to check convergence. More information can be found in Supplementary Material 2.

Species delimitation
We tested species hypotheses within Camallanidae mainly focused on the present samples. Two different species discovery methods and one species validation method were used to delimit species boundaries. Some analyses required ultrametric trees; these were generated in BEAST 2.5 as previously described and additional details are also in Supplementary Material 2.
For species discovery we used the Generalized Mixed Yule Coalescent (GMYC) method and the Poisson Tree Process (PTP) method, which do not require a priori assignments regarding putative species (Pons et al. 2006;Zhang et al. 2013). The GMYC requires an ultrametric guide tree, uses ML to delimit species and estimates a transition point before which all nodes reflects species diversification events and after which all nodes represent population-coalescent process (Pons et al. 2006;Razkin et al. 2016). GMYC tests were run in the web server (http://species.h-its.org/gmyc) under the single-threshold and multiple-threshold models, for each dataset. The PTP does not require ultrametric trees; therefore, the major consensus tree generated from BI of each dataset, were used in the analysis; this method try to identify significant changes in the rate branching of the phylogenetic tree, using number of substitutions. These analyses were run in the bPTP web server (http://species.h-its.org/ptp) using the default settings (1 × 10 5 generations, thinning set to 100, burn-in 1 × 10 4 ), which are adequate for datasets with < 50 taxa according to Zhang et al. (2013). PTP run both ML and BI to support the delimited species.
To validate the species, we used *BEAST (Heled and Drummond, 2010) implemented in BEAST 2.5 to generate unrooted (in order to improve species delimitation results; see Zhang et al. 2013) species trees based on each dataset (minus outgroup sequences). *BEAST uses BI approach to generate phylogenies. We considered samples from a same species as different populations, since they are mostly originating from different host individuals (see Table 1). Substitution models were chosen based on bModelTest, molecular clock selected based on MLE from stepping stone and path sampling, Yule process species tree priors and a constant root population size model. MCMC chains were run for 1 × 10 8 generations, sampling every 2 × 10 3 generations and ¼ of initial runs discarded as burn-in. Chains were run in duplicates and inspected with Tracer to check convergence. More information can be found in Supplementary Material 2.
We also used the Automatic Barcode Gap Discovery (ABGD) method (Puillandre et al. 2012) to generate pairwise (patristic) distance (P) matrixes, and thus evaluate intra and interspecific divergences among samples according to each genetic marker. The ABGD analyses were run online (http://wwwabi.snv.jussieu.fr/public/abgd/) using Kimura two-parameter (K2P) (Kimura, 1980) as distance metric, with other parameters set to default. We refrain from account the partitions generated by ABGD, since intraspecific divergence of 18S and COI are poorly known and editing ABGD parameters for intraspecific genetic distance would generate mistaken results.

Results
The prevalence of P. (S.) inopinatus found in the present study was 30% (3 parasitized/10 examined) in A. passionis and 25% (3/12) in M. elongatus. The mean intensity of infection was 7.3 ± 1.5 (ranging from 6-9 parasites per infected host) in A.  wider than those non gravid (both p < 0.001), the tail to body length ratio was higher in non-gravid (p < 0.001) and vulva was more posterior in gravid specimens (p = 0.01); however, the relative position of vulva was similar in both gravid and non-gravid females (p = 0.2) ( Table 2, Supplementary Material 5). Width of buccal capsule and length of oesophagus were also higher in gravid females (p = 0.014 and < 0.001, respectively); distance from anterior end to deirid and nerve ring to total body length ratio also showed statistical differences according to female ontogeny ( nodes (e.g., for comparing genera) (Fig. 3).
Regarding the species delimitation analyses, PTP was more concordant with the molecular phylogenetic approach than GMYC (Figs. 2, 3). ML of the GMYC was  (Fig. 2). In the tree of COI sequences, PTP agreed with taxonomic labelling of all species, except B. slomei and B. xenopodis that were recognised as a single species (Fig. 2), and the K2P distance between these sequences (MG948463/ MN523681) was null (P = 0). The concordance of these methods (mainly PTP) with the taxonomic labelling of the sequences was higher in the dataset of COI and 18S + COI sequences, as expected (Figs. 2, 3).
Species tree approach validated the results of phylogenetic and species delimitation methods (mainly of PTP) (Figs. 2-4). Samples of P. (S.) inopinatus were validated as a single species, whereas those that appeared as non-monophyletic were not validated as single species, and, as same as in species delimitation approach, B. slomei and B. xenopodis and recognised as a single lineage; nodal supports in these cases were medium to high (Fig. 4). We indicated host, geographic origin and environmental features related to samples in the species tree of 18S and COI (Fig. 4). Thus, most phylogenetic assemblages did not show similarities regarding these characteristics, except the clade A in Fig. 4 that grouped species from freshwater hosts and the shaded one, composed of species from India, with no additional information. These clades showed average nodal supports (Fig. 4). It should be mentioned that P.

Discussion
As previously commented, Procamallanus (S.) inopinatus is a nematode parasite with suggestively low host specificity and wide distribution in the freshwaters of South America (from North to South-western Brazil, Venezuela, Peru and Argentina) (see Moravec, 1998;Chemes and Takemoto, 2011;Luque et al. 2011Luque et al. , 2016. Some aspects of its taxonomic assessment have been controversial, especially considering the spiral ridges in buccal capsule (see Moravec et al. 1993Moravec et al. , 1997. Despite the wide variation in the number of spiral ridges in P. (S.) inopinatus (13-23, see Supplementary Material 4), they occupied no more than 2/3 of the buccal capsule as shown herein and in other studies (Travassos, 1928;Pereira, 1935;Kloss, 1966;Petter and Thatcher, 1988;Moravec et al. 1993Moravec et al. , 1997Abdalla et al. 2012). Reciprocally, Petter and Dlouhy (1985) described this structure completely occupied by the spiral ridges and Vicentin et al.
(2013) reported the same characteristic, but based only on female specimens; both studies analysed characiform fishes from River Paraguay basin (as in the present study).
These observations indicate that this characteristic may be a rare intraspecific variation in P. (S.) inopinatus. However, to confirm this assertion, specimens with the whole buccal capsule occupied by spiral ridges need to be genetically characterised and compared to the present data. It should be mentioned that, based on genetic evidence, the morphology of buccal capsule seems to bears some degree of artificiality in the systematics of genera and subgenera of Camallanidae (Wijová et al. 2006 The present results showed that almost all measurements used for taxonomic diagnosis had medium to strong correlation with the body length of parasites.
Measurements of the organs showed a positive correlation (except the length of spicules and the distance from cephalic and to vulva) whereas the ratio of these measurements to total body length showed negative correlation (see Supplementary Material 5). It indicates that the growth rates of body and organs are markedly different, in which the first grows much faster than the second (a matter of mathematic reason: if the divisor, i.e., body length, is much higher than the dividend, i.e., measure of a given organ, the quotient will be low). Moreover, these results indicate general morphometric dependence on body length, except by the spicule length and the relative position of vulva, which are more constant and, consequently, more reliable for intraspecific comparisons, regardless of body length and sexual maturation of P. (S.) inopinatus.
Sexual maturation of females of P. (S.) inopinatus in the present study had significant effect on their body length, in which gravid specimens were larger than the non-gravid. As a consequence, some ratio measurements, the maximum body and buccal capsule width, oesophagus length, distance from cephalic end to vulva, but not its relative position (vulva to body length ratio), were statistically different ( Table 2 basins (Petter and Dlouhy, 1985;Petter and Thatcher, 1988;Moravec et al. 1993Moravec et al. , 1997. highlighted that GMYC was more conservative than PTP, regarding species delimitation, and the accuracy of the methods was higher in the databases including COI (see Supplementary Material 6), because this genetic region was more informative for analysing divergences and similarities among closely related taxa than the 18S (see Fig. 3); this is due to the fact that, in metazoans, mtDNA genes evolve at faster rates than those form nuclear DNA (Allio et al. 2017).

Although the body length is similar between males
The overall phylogeny of Camallanidae was poorly resolved and monophyly was not observed in most genera and subgenera as a consequence of the morphologybased systematics that predominantly uses the buccal capsule structure as diagnostic feature, which seems to be artificial as commented previously (Wijová et al. 2006;Černotíková et al. 2011;Sardella et al. 2017;Ailán-Choke et al. 2019).
The genus Batrachocamallanus Jackson & Tinsley, 1995 is a particular case; it has been considered synonym of Procamallanus (Procamallanus) Baylis, 1923 by some authors (Moravec et al. 2006) and valid by others (Jackson and Tinsley, 1995;Svitin et al. 2018Svitin et al. , 2019. The buccal capsule structure in both taxa is the same and the characters used in the differential diagnosis of Batrachocamallanus were refuted and proved to be also present in Procamallanus (Procamallanus) (see Moravec et al. 2006). Therefore, we agree with the proposal of Moravec et al. (2006) and consider Batrachocamallanus junior synonym of Procamallanus (Procamallanus). In this sense, the present genetic analyses support the previous assertion because GMYC considered the two representatives, labelled B. slomei Southwell &Kirshner, 1937 andB. xenopodis (Baylis, 1929), along with P. (P.) spiculocubernaculus Agarwal, 1958 as a single species and in the species tree they were closely related with full nodal support.
Furthermore, sequences from B. slomei and B. xenopodis were 100% similar, but we will not assume their synonymy, since data represents only a small and partial fragment of a gene.
These relationships between habitat, host taxa, geographic origin and the phylogenetic aspects of camallanids cannot be generalised and are unclear, since as same as some assemblages were formed by samples isolated from freshwater hosts, they were from very different geographic origins and host taxa; other assemblages were formed independent from these characteristics. Similar results were previously observed     superscript numbers make correspondence with information in to the Bayesian posterior probability (BPP) and maximum likelihood (ML) bootstrap replications as follows: full squares (BPP = 1, ML = 100%), full circles (0.96 < BPP < 1, 90% < ML < 100%) and empty circles (0.90 < BPP < 0.96, 95% < ML < 90%). Black polygons next to taxa labels indicate specific entities recognised by the Poisson Tree Process (PTP) and those white by the Generalized Mixed Yule Coalescent (GMYC); superscript numbers make correspondence with information in Table 1. Sequences from de present study are in bold and superscriptions are RM = River Miranda, RX = River Xingu. follows: full squares (BPP = 1), full circles (0.96 < BPP < 1) and empty circles (0.90 < BPP < 0.96). Superscript numbers make correspondence with information in Table 1.
The shaded clade represents atypical phylogenetic behaviour and includes sequences with apparently taxonomic misidentification. Letter A highlights the clade formed exclusively by species from freshwater hosts. Order of hosts, geographic origin and environment associated with each sequence are depicted. Sequences from de present study are in bold and superscriptions are RM = River Miranda, RX = River Xingu.