Introduction
The species belonging to the genus Anisakis constitute one of the most widespread groups of ascaridoid nematodes in the marine ecosystem. These parasites involve various host species at different levels across the food webs: crustaceans as first intermediate host, fish, and squids acting as intermediate and/or paratenic ones and, finally, cetaceans as their main definitive hosts. Thus, their life-cycle completion is depending on the stability of marine trophic webs (Mattiucci and Nascetti, Reference Mattiucci and Nascetti2008; Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a).
To date, nine species belonging to the genus Anisakis have been genetically detected and morphologically described, worldwide. All those species possess distinct gene pools and are reproductively isolated; their existence as distinct phylogenetic units was demonstrated by various phylogenetic analyses (Valentini et al., Reference Valentini, Mattiucci, Bondanelli, Webb, Mignucci-Giannone, Colom-Llavina and Nascetti2006; Cavallero et al., Reference Cavallero, Nadler, Paggi, Barros and D'Amelio2011; Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014). They have been also characterized on the basis of their different ecological features, having differential geographical distribution and host ranges (for a review see Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a).
Furthermore, two more gene pools, whose nomenclatural designation and formal description still need to be defined, have been detected genetically in the genus Anisakis – provisionally indicated as Anisakis sp. 1 and Anisakis sp. 2 – (Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a).
Three closely related taxa are so far included in the A. simplex (s. l.) complex: i.e. A. pegreffii, A. simplex (s. s.) and A. berlandi (Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014). They are often distributed in metapopulations of their intermediate/paratenic (mainly pelagic/demersal fish and squids) and definitive (mainly Delphinoidea cetaceans) hosts from different geographical areas within the range of distribution of these parasites species. They have been found in co-infection in intermediate/transport hosts, such as several fish species and squids as well as in cetaceans from sympatric areas of distribution of the three sibling species (for a review, Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a). The third stage larvae of these Anisakis sibling species are also infecting a number of fish species of high commercial-economic importance (Levsen et al., Reference Levsen, Svanevik, Cipriani, Mattiucci, Gay, Hastiee, Pierce, Bušelić, Mladineo, Karl, Ostermeyer, Buchmann, Højgaard, González and Pascual2018). Humans may be accidental host, acquiring infection through consumptions of raw, marinated or undercooked fish and squids, infected by alive Anisakis spp. third stage larvae. Up to now two species, i.e. A. pegreffii and A. simplex (s. s.), are recognised as zoonotic parasites to humans (Umehara et al., Reference Umehara, Kawakami, Araki and Uchida2007; Mattiucci et al., Reference Mattiucci, Paoletti, Colantoni, Carbone, Gaeta, Proietti, Frattaroli, Fazii, Bruschi and Nascetti2017, Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a; Guardone et al., Reference Guardone, Armani, Nucera, Costanzo, Mattiucci and Bruschi2018). Worldwide, there is an increasing recognition that the fish-borne zoonosis they provoke, known as anisakiasis, is an important emerging human disease (Chai et al., Reference Chai, Murrell and Lymbery2005; Buchmann and Mehrdana, Reference Buchmann and Mehrdana2016; Guardone et al., Reference Guardone, Armani, Nucera, Costanzo, Mattiucci and Bruschi2018; Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a). Risk factors for the transmission of the zoonotic species include consumption of locally harvested wild fish and traditional seafood home preparations. There is also an increasing interest in prioritizing fish-borne zoonoses, not only in terms of their impact on human health but also in predicting the effects of climate and marine ecosystem changes. Thus, changes related to the ecology and epidemiology of actual Anisakis species may indirectly point at the consequences such changes may inflict on their natural hosts (both definitive and intermediate/paratenic ones) and the environment they live in (Mattiucci and Nascetti, Reference Mattiucci and Nascetti2008; Zarlenga et al., Reference Zarlenga, Hoberg, Rosenthal, Mattiucci and Nascetti2014; Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a).
Therefore, accurate identification of Anisakis species is essential for understanding their distribution, epidemiology, ecology, as well as their zoonotic potential in humans. It has been recognised that a multilocus genotyping approach is of great importance for species identification within the A. simplex (s. l.) complex (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016, Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a). Applying the most appropriate sets of genetic/molecular markers in Anisakis species permits to: (1) disentangle cryptic species; (2) identify them at any stage of the life cycle; (3) establish phylogenetic relationships between related taxa; (4) provide an intraspecific population analysis; (5) estimate their genetic diversity and (6) assess possible hybridization and introgression phenomena between related species. Several genetic/molecular markers have so far been tested and validated in species recognition and to study the genetic structure of Anisakis spp. (Nascetti et al., Reference Nascetti, Paggi, Orecchia, Smith, Mattiucci and Bullini1986; D'Amelio et al., Reference D'Amelio, Mathiopoulos, Santos, Pugachev, Webb, Picanḉo and Paggi2000; Valentini et al., Reference Valentini, Mattiucci, Bondanelli, Webb, Mignucci-Giannone, Colom-Llavina and Nascetti2006; Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014, Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016; Mladineo et al., Reference Mladineo, Trumbić, Radonić, Vrbatović, Hrabar and Bušelić2017). Few attempts have been made to infer information on the genetic variation at the nuclear level of A. pegreffii, A. simplex (s. s.) and A. berlandi over a large range of populations and geographical areas. Indeed, investigations of the genetic structure of the species at population level was mostly based on allozymes (Mattiucci et al., Reference Mattiucci, Nascetti, Cianchi, Paggi, Arduino, Margolis, Brattey, Webb, D'Amelio, Orecchia and Bullini1997; Mattiucci and Nascetti, Reference Mattiucci and Nascetti2007), and, more recently, partially studied by means of mitochondrial markers (Baldwin et al., Reference Baldwin, Rew, Johansson, Banks and Jacobson2011; Blažeković et al., Reference Blažeković, Pleić, Đuras, Gomerčić and Mladineo2015; Mattiucci et al., Reference Mattiucci, Giulietti, Paoletti, Cipriani, Gay, Levsen, Paoletti and Nascetti2018b). The paucity of studies in Anisakis spp. at the population level was probably due to the perceived difficulty in obtaining informative nuclear markers for the recognition of these parasites at both species and population level.
A first attempt to isolate a panel of microsatellites, also known as simple sequence repeats (SSRs), for the species A. simplex (s. s.) and A. pegreffii was performed by Mladineo et al. (Reference Mladineo, Trumbić, Radonić, Vrbatović, Hrabar and Bušelić2017). The aim of the present study was to develop a novel panel of microsatellite loci, and validate them over a large number of specimens of A. pegreffii and A. simplex (s. s.) collected from several allopatric metapopulations within their geographical distribution range in order to: (i) provide further nuclear markers that allow to distinguish the two species in a multi-marker nuclear genotyping approach; (ii) evaluate the genetic variation of these parasite species inferred from SSRs analysis; (iii) estimate the genetic differentiation between the two species of the A. simplex (s. l.) complex, i.e. A. pegreffii and A. simplex (s. s.), as inferred from the SSRs and (iv) infer data on their population genetic structure.
Materials and methods
Collection of samples
A total of N = 943 specimens belonging to the two species A. pegreffii and A. simplex (s. s) were examined at DNA microsatellite loci. The specimens were collected as L3 larval stage from intermediate/paratenic fish hosts and as L4-stage larvae and adults from cetacean species in sampling localities of the North East Atlantic Ocean, Mediterranean Sea, South West Atlantic Ocean and South West Pacific Ocean (Fig. 1 and Table 1). These localities are included in the geographical range of two species A. pegreffii and A. simplex (s. s.), from where they are reported in allopatry, according to Mattiucci et al. (Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a). Details concerning the sampling localities of the intermediate/paratenic (fish) and definitive hosts (cetaceans) of the two Anisakis species examined in this study are given in Table 1. Some of the nematodes were obtained from the frozen collection of anisakids stored in the Section of Parasitology, Department of Public Health and Infectious Diseases of ‘Sapienza University in Rome’; whereas other L3 samples were collected during the years 2013–2016 from fish in the framework of the Project PARASITE and belong to the PARASITE-Biobank (Gonzáles et al., Reference González, Gracia, Miniño, Romón, Larsson, Maroto, Regueira and Pascual2018). Node at the Section of Parasitology (Sapienza-University). The collection of Anisakis spp. adults was undertaken during the years 2011–2012 from stranded cetaceans, i.e. the pilot whale, Globicephala melas Traill, the Risso's dolphin Grampus griseus (Cuvier) and the Hector's dolphin Cephalorhynchus hectori (van Beneden) in the Southern Pacific Ocean (off the coast of New Zealand), the striped dolphin, Stenella coeruleoalba (Meyen), stranded off the Italian coast, and the harbour porpoise, Phocoena phocoena (Linnaeus) from the English Channel (Table 1).
Collecting sites of populations of A. simplex (s. s.) and A. pegreffii analysed in the present study, mapped into the geographical range of the two parasites species (for location reference codes see Table 1).

Number of specimens of A. pegreffii and A. simplex (s. s.) analysed at the eight microsatellite loci and EF1 α-1 nDNA gene, reported with host species and sampling area

N A♂♂, number of adult male nematodes; N A♀♀, number of adult female nematodes; N L4 and N L3, number of specimens at fourth and third larval stages, respectively.
Nematodes collected from the stomach of their hosts were washed in saline solution and then preserved in 70% alcohol. The nematodes collected from cetacean definitive hosts were first distinguished as L4-stage larvae and adults, then adults between females and males, according to the main morphological features diagnostic between sexes (Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a), by the use of an optical microscope at ×100–400 total magnification. The central part of the worm's body was thus used for the molecular analysis, while the cephalic and caudal ends were stored for male and female discrimination.
DNA extraction
Total DNA was extracted from each Anisakis larva L3 in toto, and from a tissue portion of 2 mg from each adult and L4 larval specimen (Table 1). The cetyltrimethylammonium bromide method (CTAB) was used as the extraction method (details in Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014). DNA obtained was quantified by using the Qubit™ dsDNA HS Assay Kit with Qubit 2.0 (Invitrogen™) (Sambrook and Russell, Reference Sambrook and Russell2001). First, aliquots with a standard concentration of 10 ng ml−1 were prepared from a pooled DNA extraction of N = 10 specimens of A. pegreffii for developing the genome library, and then on 15 unrelated specimens of the same species for testing the new markers. Later on, DNA from each specimen considered in this study was used to perform the SSRs analysis.
Genetic identification of Anisakis spp. specimens
All the specimens used (N = 943) to validate the novel DNA microsatellite loci here developed were previously identified at the species level by molecular identification, using sequences analysis at the EF1 α-1 nDNA (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016).
The amplification of EF1 α-1 nDNA was performed using the primers EF-F (5′-TCC TCA AGC GTT GTT ATC TGT T-3′) and EF-R (5′-AGT TTT GCC ACT AGC GGT TCC-3′) (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016) under the following conditions: 94 °C for 3 min (initial denaturation), followed by 35 cycles at 95 °C for 45 s (denaturation), 58 °C for 40 s (annealing) and 72 °C for 1 min, with a final post-amplification step at 72 °C for 10 min.
PCR amplification was carried out in a volume of 25 µL containing 0.5 µL of each primer 10 mm, 2.5 µL of MgCl2 25 mm (Promega), 1.5 µL of PCR buffer 5× (Promega), 0.08 mm of DMSO, 0.5 µL of dNTPs 10 mm (Promega), 5 U of HotStart Go-Taq Polymerase (Promega) and 2 µL of total DNA (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016), for the EF1 α-1 nuclear gene.
Microsatellite (SSRs markers) development and primers selection
For the development of SSRs markers, using the whole-genome sequencing via Illumina (MiSeq platform), a sub-samples of pooled 10 specimens of the species A. pegreffii was first selected for the genomic library development. In order to avoid parasite tissue contamination with the fish host DNA, alive A. pegreffii larvae collected from the fish were first washed using physiological solution, plus acetic acid at 4%. Then, the larvae were maintained in vitro culture in a Petri dish, with physiological solution for 24 h, at 6–7 °C; the larvae were thus transferred in single 1.5 mL tubes, containing PBS 1× solution, for 72 h (the PBS 1× solution was changed daily); finally, Anisakis larvae were washed several times with ultrapure water, and stored for the following DNA extraction. ‘Pure’ DNA from those A. pegreffii treated larvae was extracted and first tested by using the barcode primers for several fish species, to be sure that DNA from the parasite was not contaminated with that from the host.
Size-selected fragments from genomic DNA of the considered species (A. pegreffii) were enriched for SSRs content by using magnetic streptavidin beads and biotin-labelled GTAT, GATA, repeat oligonucleotides. The SSRs-enriched library was analysed on an Illumina MiSeq platform using the Nano 2 × 250 v2 format. After assembly, 19.332 contigs/singlets were screened and 785 contigs/singlets contained a microsatellite insert with a tetra- or a trinucleotide of at least six repeat units, were considered; finally, suitable primer design was possible in 429 microsatellite candidates. The SSR-enriched library obtained from those pooled samples of A. pegreffii was also checked with the genome of A. simplex (s. l.) (PRJEB496) hosted at WormBase Parasite (Howe et al., Reference Howe, Bolt, Cain, Chan, Chen, Davis, Done, Down, Gao, Grove, Harris, Kishore, Lee, Lomax, Li, Muller, Nakamura, Nuin, Paulini, Raciti, Schindelman, Stanley, Tuli, Van Auken, Wang, Wang, Williams, Wright, Yook, Berriman, Kersey, Schedl, Stein and Sternberg2016), assembled by the Parasite Genomic Group at the Welcome Trust Sanger Institute, UK (http://parasite.worm-base.org/Anisakis_simplex_prjeb496/Info/Index).
A total of seven A. pegreffii microsatellite loci were selected and designated as: Anisl 05784, Anisl 08059, Anisl 00875, Anisl 07132, Anisl 00314, Anisl 10535, Anisl 00185 (Table 2). Primers were designed from the flanking sequences to produce a 69–250 bp amplicons and to have an optimal annealing temperature near 56 °C (Table 2). Alleles at each locus are referred to by their size in base pairs, as determined by GeneMapper software (Table 2). The primers for the seven best loci were chosen and scored, at first, on 15 specimens of A. pegreffii collected from different fish hosts and geographical areas, for testing them in a larger set of individuals and for consistency in PCR amplification. All markers revealed a number of alleles ranging between 5 and 12 in those first analysed specimens of A. pegreffii used for SSRs development. The allele size was determined by the use of an ABI3730 instrument, using GeneScanTM-500 LIZ Size Standard.
Locus name, primers sequences, number of alleles (N A), repeat type, fragment bp size and fluorescent dye of seven microsatellite loci developed in A. pegreffii, which cross amplify in A. simplex (s. s.)

* = labelled primer.
The above-selected novel SSRs loci (Table 2) were then used to evaluate their success in cross-species amplification of A. simplex (s. s.) species. The markers developed for the specimens of A. pegreffii successfully amplified DNA obtained from A. simplex (s. s.) samples. Therefore, all the larval and adult specimens corresponding to the two species A. pegreffii and A. simplex (s. s.) reported in Table 1 were analysed with the novel panel of SSRs markers.
PCR amplification of SSRs markers, and genotyping
Two multiplex PCR amplifications (with four and three couple of primers respectively) were optimized to be performed in a 10 µl reaction volume, containing 5–10 ng of genomic DNA, 5 µl Type-it Microsatellite PCR Kit (Qiagen®), double distilled water, and concentrations of 10 µM labelled forward and reverse primers each (Table 2). The following cycling protocols were used for the amplification both for the two multiplex reactions: 35 cycles with 94 °C for 30 s, 56 °C for 90 s, and 72 °C for 60 s. Before the first cycle, a prolonged denaturation step (95 °C for 15 min) was included and the last cycle was followed by a 15 min extension at 60 °C. The features for the seven loci selected, i.e. the repeat motif, the bp size, the primer sequences and its labelling dye, are reported in Table 2.
Amplified PCR products were genotyped by an external Company (Macrogen service).
An additional microsatellite marker, i.e. the Anisl 7 from Mladineo et al. (Reference Mladineo, Trumbić, Radonić, Vrbatović, Hrabar and Bušelić2017), was also tested on all those specimens of the two Anisakis spp. considered in the present study (Table 1); its analysis was included in the Multiplex no. 2 and the forward primer was labelled with the NED dye (Table 1); specific PCR conditions used for this locus were those previously described in the text for multiplex 2 (see above).
Finally, the alleles obtained from the electropherograms were identified and binned using the software Genemapper v.4.1 (Applied Biosystems, USA). Genotyping errors generally associated with microsatellite analysis, such as stutter bands, the presence of null alleles and large drop-out alleles, were verified by using the software MICRO-CHECKER version 2.2.3 (Van Oosterhout et al., Reference Van Oosterhout, Hutchinson, Wills and Shipley2004). Patterns of tri- and tetrallelic peaks in the female individuals, as possible results of tissue contamination from sperm males fecundation, were not found.
Genetic data analysis
All samples tested in this study at the panel of microsatellite markers were previously identified by sequences analysis obtained at the EF1 α-1 nDNA partial gene (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016). The sequences obtained in the EF1 α-1 nDNA region were aligned by using Clustal X (version 2.0) software (Larkin et al., Reference Larkin, Blackshields, Brown, Chenna, McGettigan, McWilliam, Valentin, Wallace, Wilm, Lopez, Thompson, Gibson and Higgins2007) with those from the same species previously obtained (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016), in order to detect those fixed diagnostic nucleotide positions allowing to discriminate A. simplex (s. s.) and A. pegreffii (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016).
The observed heterozygosity (H o), the expected heterozygosity (H e), the mean number of alleles per locus (A) were calculated by ARLEQUIN version 3.5 (Excoffier and Lischer, Reference Excoffier and Lischer2010). The occurrence of the expected Hardy–Weinberg equilibrium (HWE) at each locus was assessed by means of exact tests, as implemented in the software ARLEQUIN version 3.5 (Excoffier and Lischer, Reference Excoffier and Lischer2010). Significance levels were adjusted using the sequential Bonferroni correction for multiple tests (Rice, Reference Rice1989). The fixation index (F IS) from the genetic data sets obtained at the SSRs loci in the analysed populations of the two Anisakis species, was estimated using ARLEQUIN version 3.5 (Excoffier and Lischer, Reference Excoffier and Lischer2010), and F st (Weir and Cockerham, Reference Weir and Cockerham1984) by FSTAT version 1.2 (Goudet, Reference Goudet1995). We used analysis of molecular variance (AMOVA) as implemented in ARLEQUIN version 3.5, with 1000 permutations, on the genetic data sets obtained in the populations of A. simplex (s. s.) and A. pegreffii. Nei's distance matrix (Nei, Reference Nei1978), was calculated from the allele frequencies estimates, by using BIOSYS 2.0 software program (Swofford and Selander, Reference Swofford and Selander1997). An UPGMA (Unweighted Pair Group Method with Arithmetic Mean) was generated by using PHYLIP software (Felsenstein, Reference Felsenstein1993), based on Nei's distance values. In addition, a multivariate-based method for inference of population structuring was based on non-metric multidimensional scaling (nMDS) ordination plot, by using allele frequencies observed at the microsatellite loci in all the tested populations of the two parasites species, based on the calculation of Euclidean distance. The nMDS-plot for the data set was performed by the package Vegan Community Ecology Package, R version 2.5-4 (Oksanen, et al., Reference Oksanen, Blanchet, Friendly, Kindt, Legendre, McGlinn, Minchin, O'Hara, Simpson, Solymos, Stevens, Szoecs and Wagner2019).
Finally, a Bayesian clustering algorithm elaborated by the program STRUCTURE 2·3·3 (Pritchard et al., Reference Pritchard, Stephens and Donnelly2000) was used to assign each analysed specimen to the species A. pegreffii or A. simplex (s. s.), as based on its genotype showed at those SSRs markers, plus the EF1 α-1 nDNA region. The same approach allowed to identify any eventual instances of gene exchange between the two species. STRUCTURE is a model-based procedure that uses individual multi-locus genotypes to identify the optimal number of clusters (K) in a dataset, by minimizing the resulting Hardy-Weinberg and linkage disequilibria. The analysis was run setting the predefined number of clusters between 1 and 9 (i.e. the number of sampling areas in our dataset). Twenty replicates of the analysis were carried out to check for consistency; each run for 100 000 MCMC iterations, following a burn-in of 50 000 iterations, under the admixture model and the assumption of correlated allele frequencies among populations. The best K value was identified using both the log probability of the data and the rate of change in the log probability of the data between successive K values as optimality criteria (Evanno et al., Reference Evanno, Regnaut and Goudet2005).
Results
Identification of Anisakis spp. specimens by sequencing of the EF1 α-1 nDNA region
A total of N = 943 larval and adult specimens of Anisakis spp. (Table 1) to be used for the microsatellites development and validation, were first identified to species level by sequences analysis of the nuclear EF1 α-1 nDNA partial region. A fragment of 409 bp of the EF1 α-1 nDNA region was obtained for all 943 specimens analysed. It revealed the presence of two diagnostic nucleotide sites in positions were 186 and 286 (as indicated in Mattiucci et al. (Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016)), showing T and C in A. pegreffii but C and T in the parental taxon A. simplex (s. s.), respectively. According to this analysis, N = 451 specimens were assigned to species A. pegreffii, whereas N = 492 were assigned to A. simplex (s. s.) (Table 1). Patterns of heterozygote genotypes, i.e. showing two overlapping C/T peaks, were not found among all the presently analysed individuals of the two species. Sequences of the EF1 α-1 nDNA region were deposited in GenBank under the accession numbers from MK032267 to MK032270 for A. pegreffii, and from MK032271 to MK032275, for A. simplex (s. s.).
Genetic diversity of SSRs in A. pegreffii and A. simplex (s. s.)
The number of individual specimens genotyped at the eight scored SSRs varied from a minimum of N = 77 to a maximum N = 164 worms of A. pegreffii and from N = 57 to N = 176 from the selected populations of A. simplex (s. s.) (Table 3). For some markers, such as the locus Anisl 7, there was a small number of worms in each population of the two analysed parasite species that repeatedly failed in the genotyping. Thus, the total number of nematodes of the two species tested at the locus Anisl 7, was N = 713, with respect to the total sample (N = 943) (Tables 1 and 3). However, a high representative number of specimens from each population of the two species from the selected geographic areas was available at that locus to be used in the subsequent genetic data analyses (Table 3).
Genetic diversity at eight microsatellite loci of adult and larval specimens in nine populations of A. pegreffii and A. simplex (s. s.) analysed in the present study

For population codes see Table 1. N, total number of genotyped nematodes at each locus; H o, observed heterozygosity; H e, expected heterozygosity; A, number of alleles detected at each locus; U, unique alleles, occurring at frequency >0.1, among the populations of the two species. P, indicates the significance (P < 0.05) of the deviation from HWE expectation. ***P ≪ 0.001, *P < 0.01.
There were some differences in the overall genetic diversity values of the different populations of the two species, based on both the mean number of alleles per locus, and values of expected heterozygosity per locus (H e) (Table 3). All the microsatellite loci scored in the present study were polymorphic, with the total number of alleles varying in A. pegreffii between 2 (such as those observed at the locus Anisl 7 in the population from the Adriatic Sea (AD)) and 13 (such as those found in the population from New Zealand (NZ) at the locus Anisl 05784) (Table 3). In A. pegreffii, the mean value of alleles per locus (A) was A = 8.21. In A. simplex (s. s.), the total number of alleles was found to vary between 3 (such as those found at the locus Anisl 10535 in the population from Baltic Sea, BA) and 13 (such as those observed at the loci Anisl 05784 and Anisl 00314 in the population from English Channel, EC) (Table 3). In A. simplex (s. s.), the mean value of alleles per locus (A) was A = 8.27.
Deviations from the Hardy-Weinberg Equilibrium (HWE) at each locus were tested in each selected population of the two species (Table 3). Generally, positive values of F IS indicated an excess of homozygote genotypes at the selected loci; whereas, negative values indicate an excess of heterozygotes from the expected HWE (Fig. 2A). No significant deviations between observed (H o) and expected (H e) heterozygosity resulted at loci Anisl 10535, Anisl 07132, Anisl 05784, Anisl 08059 and Anisl 00875, in all the analysed populations of A. pegreffii and A. simplex (s. s), with not-significant departure from the HWE (Table 3). Some significant HWE deviations occurred in A. simplex (s. s.) at the locus Anisl 00185; conversely, the same locus did not show significant deviation between H o and H e in the examined populations of A. pegreffii (Table 3) (Fig. 2A). However, due to the high significant departure from HWE observed at that locus in populations of the species A. simplex (s. s.), this marker was not taken into account when population genetic data have been considered.
(A) F IS calculated at the eight microsatellites loci studied in the two species A. pegreffii and A. simplex (s. s.). Negative values indicate heterozygous excess while positive values indicate homozygous excess from that expected under Hardy–Weinberg Equilibrium (HWE) and (B) F IS in male and female specimens of A. pegreffii and A. simplex (s. s.), at the two sex-linked loci, i.e. Anisl 00314 and Anisl 7.

Sex-linked microsatellite loci
Statistically high significant departures from the HWE in all the analysed populations of A. pegreffii and A. simplex (s. s.) were observed also at loci Anisl 00314 and Anisl 7. The two loci showed a marked excess of homozygotes in all populations (Table 3), with positive F IS value (Fig. 2A). Interestingly, when splitting the genotypes at the loci Anisl 00314 and Anisl 7 found in both A. pegreffii and A. simplex (s. s.) into adult male and female worms, it resulted that the male worms of the two species were homozygous at those loci (F IS = 1) (Fig. 2B). Therefore, the loci Anisl 00314 and Anisl 7 appear to be sex-linked loci in both species, likely because of the hemizygosity state (i.e. diploid individuals in which there is only one allele present at the sex-linked loci) of the males in both the Anisakis spp., as X-linked inheritance. Indeed, no significantly departures from the HWE were found between the expected homozygous and heterozygous females at the Anisl 7 in the two Anisakis spp. (Table 4). Analogously, no significant departure from the HWE was found at the locus Anisl 00314 in all the populations of A. pegreffii (Table 4). While, at the same locus, the populations of A. simplex (s. s.) showed a significant excess of homozygous genotypes, also after accounting only the genotypes from female worms (Table 4). The reason for this partial departure is not known, but a possibility could be that some alleles failed to amplify in the species A. simplex (s. s.), thereby increasing the apparent homozygosity of the female populations, or probably the phenomenon was due to null alleles in A. simplex (s. s.) at this locus. As a consequence of the HWE departures found in A. simplex (s. s.) (Table 4), the locus Anisl 00314 was excluded from the subsequent population genetic analysis.
Observed (H o) and expected (H e) heterozygosity at the two sex-linked loci (i.e. Anisl 00314 and Anisl 7) in A. pegreffii and A. simplex (s. s.), estimated in adult female specimens

N = number of individuals genotyped at those loci. For the population codes, see Table 1. P = indicates the significance (P < 0.05) of the deviation from HWE expectation. *P < 0.01, **P < 0.0001, ***P < <0.0001.
Whereas, to not exclude genetic data set obtained at the other sex-linked locus (i.e. Anisl 7) from the subsequent population genetic analysis, the most reliable estimates of allele frequencies of the two parasite species were calculated only in adult specimens, according to the sex-linked genetic model estimate, assuming: (i) the hemizygosity of the males in the two Anisakis species; (ii) their adult female counterparts, as biallelic nematodes at the sex-linked loci. Thus, only adult nematodes of A. pegreffii collected in S. coeruleoalba from Adriatic Sea (AD), in G. melas and C. hectori from New Zealand waters, as well as adult worms of A. simplex (s. s.) recovered in Phocoena phocaena from English Channel (EC), were accounted for the allele frequencies estimates (Table 5). Indeed, because the gender assignment was not possible for larval individuals, as a consequence, the populations of the two Anisakis species including only larvae (see Table 1) were not included in the allele frequencies estimate.
Allele frequencies observed at six microsatellite loci tested in the nine populations of A. pegreffii and A. simplex (s. s.)

For population codes, see Table 1. With regard to the sex-linked locus (i.e. Anisl 7) the most reliable estimate of allele frequencies was calculated only on adult populations of the two species, according to the sex-linked genetic model estimate, assuming: (i) the hemizygosity of the males at that locus in the two Anisakis species; (ii) their adult female counterparts, as biallelic at the sex-linked loci. Populations of the two species including only larval specimens (see Table 1) tested at the sex-linked locus Anisl 7, were excluded from the allele frequencies estimate.
Identification of Anisakis spp. specimens by using SSRs loci
Interestingly, among the SSRs loci studied, Anisl 10535, Anisl 05784, Anisl 00875, Anisl 08059 and Anisl 7 have shown a certain differential frequency distribution at some alleles in the two species (Table 5). In particular, the locus Anisl 7 had up to eight different alleles in populations of A. simplex (s. s.); it resulted to be less polymorphic in A. pegreffii, showing up to four alleles in the population from New Zealand (Table 5). However, at this locus all the alleles found in A. simplex (s. s.) populations were clearly distinct from those observed in A. pegreffii. In other words, Anisl 7 showed alternative alleles in A. pegreffii and A. simplex (s. s.), thus suggesting that this is a diagnostic locus (100%) between the two parasite species (Table 5).
Individual Anisakis worms were assigned to either A. pegreffii or A. simplex (s. s.) based on a Bayesian clustering algorithm which is implemented in STRUCTURE 2.3.3 (Pritchard et al., Reference Pritchard, Stephens and Donnelly2000), according to their genotypes observed at the five SSRs loci, i.e. Anisl 10535, Anisl 05784, Anisl 00875, Anisl 08059 and Anisl 07132 (Fig. 3A). Because of the hemizygosity of the males at the sex-linked locus Anisl 7, and the presence of null alleles found at the locus Anisl 00185 in A. simplex (s. s.), genotypes found at those loci were obviously excluded from the STRUCTURE analysis. Using both the highest ln probability and the delta-K (Evanno et al., Reference Evanno, Regnaut and Goudet2005) optimality criteria, the STRUCTURE analysis indicated K = 2 as the clustering option which best fitted the data set (Suppl. Fig. 1A and B).
Percentage contribution (Q value) of A. pegreffii and A. simplex (s. s.) to the multi-locus genotype of each studied individual (barplot), estimated at five SSRs loci (i.e. Anisl 10535, Anisl 05784, Anisl 00875, Anisl 07132 and Anisl 08059), among those here studied (A); and at the same five loci plus EF1 α-1 (B), by using STRUCTURE, with k = 2. Codes of the sampling areas (see Fig. 1 and Table 1) are as follows: NW: Norwegian Sea; NS: North Sea; GS: Grand Sole bank; EC: English Channel; BA: Baltic Sea; AD: Adriatic Sea; TY: Tyrrhenian Sea; NZ: New Zealand coast and AR: Argentine coast.

Thus, considering those five microsatellite loci with that clustering option, all individuals from samples TY (Tyrrhenian Sea), AD (Western Adriatic Sea), NZ (New Zealand coast) and AR (Argentine coast) were assigned with high percentage of assignment to A. pegreffii (Fig. 3A). Similarly, all individuals from sampling locality NS (North Sea), Norwegian Sea (NW), English Channel (EC), Baltic Sea (BA) and Grand Sole Bank (GS) were referred to A. simplex (s. s.), as inferred from the same Bayesian assignment/cluster analysis (Fig. 3A). Finally, when taking into account the genotypes of all the examined individuals found at those five nuclear loci, plus those observed at the diagnostic positions of the EF1 α-1 nDNA partial region, as clustering option by STRUCTURE, all the individuals were assigned to the distinct species (A. pegreffii or A. simplex (s. s.)), with high membership level (>99%) (Fig. 3B).
However, higher percentage of assignment (>99%) was even obtained when based on the genotypes observed only in the female specimens tested at six SSRs loci (including also the females genotyped at the X-linked locus Anisl 7, but excluding the locus Anisl 00185 due to null alleles) (Suppl. Fig. 2).
No genotypes/individuals of mixed ancestry between the two species (such as F1 hybrids), i.e. with expected Q-value = 0.50, were found in this analysis, thus suggesting that all individuals Anisakis presently analysed from allopatric populations of the two species, were “pure parental” specimens of A. pegreffii or A. simplex (s. s.). On the other hand, no heterozygotes at the diagnostic positions of EF1 α-1 nDNA, nor at the diagnostic alleles observed at the Anisl 7 locus, were found in the samples.
Genetic differentiation at the inter and intra-specific level, as inferred from SSRs loci
Excluding the SSR locus (i.e. Anisl 00185) affected by null alleles in A. simplex (s. s.), and the sex-linked loci Anisl 7 and Anisl 00314 due to hemizygosity of males in both species of Anisakis and the presence of null alleles in the latter SSR locus, the remaining loci showed adequate genetic diversity for population-level genetic analysis. The AMOVA analysis of the five nuclear markers showed that in both the species most of the variance was significantly allocated within individuals (≈87% and ≈89%), with F IT = 0.12 and F IT = 0.10, respectively, in A. pegreffii and A. simplex (s. s.). Further, a notable variation was found among individuals within populations. AMOVA provided also moderate significant variation among infrapopulations of A. pegreffii; less significant was the percentage of variation among the populations of the species A. simplex (s. s.) (Suppl. Table 1a, b).
Pairwise F st values were calculated for all the population pairs of the two Anisakis species, to quantify the extent of genetic sub-structuring between different populations (Table 6). Based on those five SSRs loci considered, a significantly high level of genetic differentiation resulted at the interspecific level, between the populations of A. pegreffii and the populations of A. simplex (s. s.) (on average F st ≈ 0.29 P ≪ 0.001) (Table 6), with the Adriatic Sea (AD) population of A. pegreffii diverging the most from the present A. simplex (s. s.) populations.
Population pairwise F st (below the diagonal) estimates by FSTAT (Goudet, Reference Goudet1995) and Nei's genetic distance (Nei, Reference Nei1978) values (D Nei, above the diagonal) calculated by BIOSYS.2 (Swofford and Selander, Reference Swofford and Selander1997) from the allele frequencies observed at five SSRs loci (excluding Anisl 00185, Anisl 00314 and Anisl 7), among the nine analysed populations of A. pegreffii (see Table 1, Fig. 1) (NZ, New Zealand coast; AR, Argentine coast; TY, Tyrrhenian Sea; AD, Adriatic Sea), and A. simplex (s. s.) (GS, Grand Sole Bank; EC, English Channel; NS, North Sea; NW, Norwegian Sea; BA, Baltic Sea)

Significant genetic differentiation: *P < 0.01; **P < 0.001), based on 1000 permutations.
At the intraspecific level, i.e. when analysing the microsatellite data separately for each species including all stages, much lower values of genetic differentiation was observed. No significant genetic differentiation was observed between conspecific populations of larvae vs adults from fish and cetacean hosts being, F st ≈ 0.008 and F st ≈ 0.004 respectively, in A. pegreffii and A. simplex (s. s.). At the interpopulation level, the genetic differentiation among populations of A. simplex (s. s.) from NE Atlantic waters was, on average, F st ≈ 0.008. The highest F st value = 0.021 (P < 0.001) was found between the Baltic Sea (BA) population of A. simplex (s. s.) vs the North Sea (NS) population of the same species (Table 6).
In A. pegreffii, low level of differentiation was found between pairs of populations geographically close to each other, such as the two samples from the Mediterranean Sea (i.e. TY vs. AD, F st = 0.002), or between the austral populations of A. pegreffii from off the New Zealand and Argentine coast (F st = 0.007) (Table 6). Whereas, higher significant genetic sub-structuring (on average, F st ≈ 0.060 P < 0.001) seems to exist when comparing the Austral populations (NZ and AR) of A. pegreffii with the Mediterranean ones (TY and AD) (Table 6).
The genetic distance (D Nei) values were based on the allele frequencies observed at five SSRs loci (Anisl 10535, Anisl 05784, Anisl 00875, Anisl 08059 and Anisl 07132). Indeed, in the estimation of the genetic distance values at population level, it was not possible to include the sex-linked locus Anisl 7 due to the fact that allele frequencies of adult specimens were available only in three populations of the two Anisakis species (i.e. NZ, AD, and EC) (see Tables 1 and 5).
At the interspecific level, the value of D Nei observed among the currently analysed populations of the two species was on average D Nei ≈ 0.91, as based on five SSRs loci. Genetic distance values between the two Anisakis spp. ranged from D Nei = 0.670 (BA vs NZ) to D Nei = 1.094, observed between the population of A. simplex (s. s.) from Norwegian Sea (NW) vs A. pegreffii from Adriatic Sea (AD) (Table 6). At the intraspecific level, values of D Nei were clearly lower. In A. pegreffii, on average, the value resulted D Nei ≈ 0.15 between the Austral and Boreal populations of the species (Table 6). Conversely, low level of genetic differentiation was found among the Mediterranean populations of A. pegreffii (D Nei = 0.005). Analogously, a high genetic similarity was observed among the NE Atlantic populations of the species A. simplex (s. s.) (on average, D Nei = 0.021), with the population from the English Channel being the most distant population (on averages D Nei = 0.037) from the other conspecific ones of the NE Atlantic Ocean (Table 6).
The genetic relationships among the populations and species of the two Anisakis taxa are illustrated in the UPGMA cluster (Fig. 4), as inferred from Nei's genetic distance values (Table 6), and in the nMDS analysis (Fig. 5) obtained from the allele frequencies calculated at five microsatellite loci (Table 5). Both analyses agreed in depicting a high level of genetic divergence between the analysed populations of A. pegreffii and A. simplex (s. s.). On the other hand, they were congruent in showing that the populations of A. pegreffii and A. simplex (s. s.) are clustering in two well distinct clades, corresponding to the two biological species. Additionally, they were also concordant in showing close genetic relationship and similarity between the two Austral populations of A. pegreffii (NZ and AR) which are clustering in the same sub-clade (UPGMA analysis) and represent close clusters in the nMDS analysis, as well. However, marked sub-structuring of the Austral populations (NZ and AR) of A. pegreffii compared to populations of the same species from Boreal areas (i.e. from Tyrrhenian Sea, TY and Western Adriatic Sea, AD) was shown in the multilocus plot analysis (Figs 4 and 5).
The unweighted pair group method of analysis (UPGMA) cluster based on Nei's genetic distance values inferred from allelic frequencies calculated at five microsatellite loci (i.e. Anisl 10535, Anisl 05784, Anisl 00875, Anisl 07132, Anisl 08059), showing the genetic relationship between populations and species here studied (. A. pegreffii and A. simplex (s.s.)). Populations and their reference codes are reported in Table 1 and Fig. 1. They are: NW: Norwegian Sea; NS: North Sea; GS: Grand Sole bank; EC: English Channel; BA: Baltic Sea; AD: Adriatic Sea; TY: Tyrrhenian Sea; NZ: New Zealand coast and AR: Argentine coast.

Non-metric multidimensional scaling (nMDS) ordination plot inferred from Euclidean distance values obtained from the allele frequencies observed at five microsatellite loci (i.e. Anisl 10535, Anisl 05784, Anisl 00875, Anisl 07132, Anisl 08059), showing the genetic relationships of populations and species here studied (i.e. A. pegreffii and A. simplex (s.s.)). The nMDS-plot for the data set was performed by Vegan Package. R package version 2.5-4 (Oksanen et al., Reference Oksanen, Blanchet, Friendly, Kindt, Legendre, McGlinn, Minchin, O'Hara, Simpson, Solymos, Stevens, Szoecs and Wagner2019). The variance percentage explained by the two axes is the following: 47.7% for axis_1 and 7.1% for axis_2; stress value = 0.0003. Populations and their reference codes are reported in Table 1 and Fig. 1. They are: EC: English Channel; GS: Grand Sole bank; NW: Norwegian Sea; NS: North Sea; BA: Baltic Sea; AD: Adriatic Sea; TY: Tyrrhenian Sea; NZ: New Zealand coast and AR: Argentine coast.

Discussion
Genetic diversity of microsatellites in A. pegreffii and A. simplex (s. s.)
Although microsatellites have been extensively used to determine population structure of marine species, so far they have received little attention in marine parasites. In recent years, microsatellites analysis of trematode parasites has been used to identify origin of the steelhead trout Oncorhynchus mykiss (Criscione et al., Reference Criscione, Cooper and Blouin2006), or to show congruent population genetic pattern of parasite species (i.e. Plagioporus shawi) and its salmonid hosts (Oncorhynchus spp.) (Criscione and Blouin, Reference Criscione and Blouin2007). Despite their potential utility, microsatellites as genetic markers for parasitic nematodes have been not much explored (Johnson et al., Reference Johnson, Webster, Adam, Buckland, Dawson and Keller2006; Criscione et al., Reference Criscione, Anderson, Raby, Sudimack, Subedi, Rai, Upadhayay, Jha, Williams-Blangero and Anderson2007; Redman et al., Reference Redman, Packard, Grillo, Smith, Jackson and Gilleard2008; Temperley et al., Reference Temperley, Webster, Adam, Keller and Johnson2009; Betson et al., Reference Betson, Halstead, Nejsum, Imison, Khamis, Sousa-Figueiredo, Rollinson and Stothard2011; Reid et al., Reference Reid, Hoareau and Bloomer2012; Patrelle et al., Reference Patrelle, Jouet, Lehrter and Ferté2014; Rabelo et al., Reference Rabelo, de Miranda, Furtado, Redondo, Tennessen and Blouin2017; Greeff et al., Reference Greeff, Reid, Gagjee, Clift and de Waal2018). This may be explained by some complexity when using microsatellites, which could arise from potential homoplasy, i.e. the independent mutation of microsatellite markers of the same size, and the occurrence of null alleles (Grillo et al., Reference Grillo, Jackson and Gilleard2006; Criscione et al., Reference Criscione, Anderson, Raby, Sudimack, Subedi, Rai, Upadhayay, Jha, Williams-Blangero and Anderson2007; Redman et al., Reference Redman, Packard, Grillo, Smith, Jackson and Gilleard2008; Glenn et al., Reference Glenn, Lance, McKee, Webster, Emery, Zerlotini, Oliveira, Rollinson and Faircloth2013).
In the present study, seven novel nuclear polymorphic loci have been developed and characterized as nuclear markers by next-generation sequencing approach in the marine parasite species A. pegreffii. The primer pairs have allowed the robust amplification from each single individual worm, and a clear unambiguous genotype detection. Further, they have cross-species amplified the same SSRs in the sister taxon A. simplex (s. s), as well as in A. berlandi (Mattiucci et al., pers.com.). However, among the novel microsatellite loci here scored, the loci Anisl 00185 and Anisl 00314 were found to be affected by null alleles in the species A. simplex (s. s.). The high frequency of null alleles generally observed in SSRs markers is likely to reflect polymorphism of flanking sequence in form of SNPs and/or indels, such as those frequently found in other nematode species (Temperley et al., Reference Temperley, Webster, Adam, Keller and Johnson2009). Also, when a mutation occurs in the primer region, some individuals could fail in one allele amplified (Selkoe and Toonen, Reference Selkoe and Toonen2006). This would have been occurred just in the species A. simplex (s. s.) because of the fact that the SSRs primers have been here first selected in the species A. pegreffii. In many other taxonomic groups, PCR amplifications are known to fail with higher likelihood when heterospecific SSRs primers are employed (Lowe et al., Reference Lowe, Jones, Raybould, Trick, Moule and Edwards2002). However, the high significant HWE deviations found at the locus Anisl 00185 and Anisl 00314 in the species A. simplex (s. s.) (Table 3 and 4) was taken into account, excluding those loci when population genetic data sets of that species were being analysed.
The novel markers detected in the present study have been validated on a very large number of individuals from allopatric populations of the two species (i.e. N = 451 of A. pegreffii and N = 492 of A. simplex (s. s.)). Further, a large number among those nematodes of A. pegreffii (N = 372) and of A. simplex (s. s.) (N = 341) were characterized at Anisl 7 SSR marker. The last one was firstly scored by Mladineo et al. (Reference Mladineo, Trumbić, Radonić, Vrbatović, Hrabar and Bušelić2017), in the aim of developing microsatellites loci in species of the A. simplex (s. l.) complex. However, the authors scored at that locus only very few specimens for A. simplex (s. s.), and from a single population of A. pegreffii from the Eastern Adriatic Sea. Thus, on the basis of results obtained on a too low number of specimens analysed by the authors, the locus Anisl 7 was not pointed out as a valuable diagnostic marker between the two species.
X-associated SSRs loci in Anisakis spp.
An interesting result was that in the species A. pegreffii and A. simplex (s. s.) the two loci Anisl 00314 and Anisl 7 are located on sexual X chromosome(s), thus being sex-linked. Indeed, it was demonstrated that male specimens of A. pegreffii and A. simplex (s. s.) are hemizygous at both the SSRs loci for several distinct alleles. However, it was also found that the locus Anisl 00314 was affected by null alleles in A. simplex (s. s.), because of the significant HWE departures observed, even when excluding male genotypes from the test (Table 4).
The finding of X-associated SSRs loci was not detected by Mladineo et al. (Reference Mladineo, Trumbić, Radonić, Vrbatović, Hrabar and Bušelić2017), who retained that the occurrence of excess of homozygote genotypes found at the same Anisl 7 locus in the two Anisakis species was related to the fact that it was affected by null alleles.
The finding of two sex-linked nuclear SSRs loci in species of the genus Anisakis represents the first evidence in anisakid nematodes. The X-linked SSRs loci have been previously reported in other ascaridoid nematodes, such as in Ascaris spp. (Criscione et al., Reference Criscione, Anderson, Raby, Sudimack, Subedi, Rai, Upadhayay, Jha, Williams-Blangero and Anderson2007), as well as several other nematodes (Johnson et al., Reference Johnson, Webster, Adam, Buckland, Dawson and Keller2006). It was indeed demonstrated that A. suum has 19 autosomes and 5X chromosomes (2n = 38A + 10X in females, 38A + 5X in males) (Müller and Tobler, Reference Müller and Tobler2000). The study on chromosomes of anisakid worms has never been carried out. However, the finding of the sex-linked loci could be due to the fact that male specimens belonging to Anisakis spp. likely possess XO sexual karyotype, as other ascarids (Müller and Tobler, Reference Müller and Tobler2000). According to Criscione et al. (Reference Criscione, Anderson, Raby, Sudimack, Subedi, Rai, Upadhayay, Jha, Williams-Blangero and Anderson2007), the number of sex-linked SSRs in the terrestrial ascarid A. lumbricoides may represent 20% of the parasite's total genome. Similarly, Johnson et al. (Reference Johnson, Webster, Adam, Buckland, Dawson and Keller2006) found that 3 of 21 SSRs loci were sex-linked in Trichostrongylus tenuis. In the present study, it seems that the percentage of SSRs X-linked loci would be around 25% in species of the A. simplex (s. l.).
Paradoxically, in spite of the hemizygosity of male nematodes, those gender-associated SSRs markers discovered in the present study could help in revealing an alternative mode of sex determination in those Anisakis species, when other gender features are not otherwise evident. Indeed, the two X-linked loci (i.e. Anisl 7 and Anisl 00314) can be considered as possible markers for the identification of an individual's female sex, at any life-history stage (even as L3 larva or immature worm), when it shows a heterozygote genotype at one or both the two loci. According to the two sex-linked loci so far discovered, sex determination would be inferred at least in the species A. pegreffii, not affected by null alleles at the locus Anisl 00314 (Table 4). However, despite their utility in gender markers, only two loci are not enough in female recognition, taking into account the loss of homozygote genotypes in females at those X-linked loci. Further X-linked SSRs loci which will be discovered in future analysis in the members of the A. simplex (s. l.) would help to increase the probability of a correct gender assignment in these nematodes.
Utility of microsatellite markers in the identification of A. pegreffii and A. simplex (s. s.)
Another aspect highlighted by this study was the finding that the present SSRs loci are of potential value in the discrimination of the two Anisakis species, when using a Bayesian clustering approach. Interestingly, in the present Bayesian analysis by STRUCTURE of genotypes observed at five loci, all the considered Anisakis specimens were assigned to A. pegreffii or A. simplex (s. s.) at very high assignment rate (>99%) (Fig. 3A). In addition, when including in the same analysis the genotypes inferred from another diagnostic nuclear marker, i.e. from the EF1 α-1 nDNA, percentage of assignment reached 100% in almost all individuals, while in others the percentage of correct assignment was, again, >99% (Fig. 3B).
Among the scored SSRs markers, it was found that some of the polymorphic loci seem to share the same alleles at several loci in A. pegreffii and A. simplex (s. s.), but often with significant different frequency proportion (such as Anisl 10535); while other loci show alternative alleles in the two species, such as the locus Anisl 7 (Table 5). Thus, the last one can be retained as a 100% diagnostic locus to separate the two species, as inferred on the present study, validated by results obtained on a large number of individuals (N = 713) of both species. A molecular/genetic key of potential value in the identification of the specimens belonging to the two species, based on the alleles scored at some SSRs loci (which can be combinable in a single one multiplex PCR analysis), has been also here proposed (Table 7).
Example of a molecular key for the identification of specimens of A. pegreffii and A. simplex (s. s.), based on the alleles detected at four SSRs loci (i.e. Anisl 05784, Anisl 08059, Anisl 10535, Anisl 7) (see Table 5), to be combinable in a single one multiplex PCR reaction

*indicates alleles found in larval samples of A. simplex (s. s.); they are not reported in Table 5, because the estimation of the allele frequencies was calculated at the sex-linked locus Anisl 7 only on adult (female and male) specimens.
In the present genetic analysis, no individuals showing evidence of mixed ancestry genotypes were detected between the two species collected from allopatric areas by the Bayesian multilocus approach. On the contrary, it has been demonstrated that misidentifications of specimens of parental A. pegreffii and A. simplex (s. s.), and their hybrid categories could result from the identification derived from only a single nuclear marker, i.e. the ITS region of the rDNA gene locus (Mattiucci et al., Reference Mattiucci, Acerra, Paoletti, Cipriani, Levsen, Webb, Canestrelli and Nascetti2016). These findings raise the question in the reliability of the ITS nuclear marker in the species identification of these two sibling species, and underline the usefulness and validity of a multilocus genotyping approach. Combining the use of several nuclear markers, such as those presently outlined, would greatly increase the informative value of future analyses in the aim to clarify patterns of hybridization and introgression events in closely related taxa of the A. simplex (s. l.) complex (Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a).
SSRs markers for population genetics analysis of A. pegreffii and A. simplex (s. s.)
The panel of microsatellite markers here identified showed to be suitable for population genetics analysis of A. simplex (s. s.) and A. pegreffii, as observed by examining the degree of genetic variation between populations of the two species collected from allopatric areas of the parasites.
Using different clustering approaches (STRUCTURE, UPGMA and nMDS analyses) as inferred from the analysis of microsatellite loci, it was demonstrated that the populations analysed in the present study, despite wide host range and geographically distant sampling localities, were strongly associated with two distinct panmictic units, that correspond to the two species A. pegreffii and A. simplex (s. s.).
At the intraspecific level, a moderate genetic differentiation (on average, F st = 0.060) was found among the Boreal populations of A. pegreffii from the Adriatic and Tyrrhenian populations of the Mediterranean Sea, and their conspecific ones from the Austral region off New Zealand and the Argentine coast. This finding is in accordance with similar values of genetic differentiation previously found at the allozyme level: indeed, in A. pegreffii a G st ≈ 0.045 value was recorded between the Adriatic Sea and New Zealand samples (Mattiucci et al., Reference Mattiucci, Nascetti, Cianchi, Paggi, Arduino, Margolis, Brattey, Webb, D'Amelio, Orecchia and Bullini1997). Also, a similar trend was observed at the mitochondrial gene (mtDNA cox-2 sequences data analysis) level (Blažeković et al., Reference Blažeković, Pleić, Đuras, Gomerčić and Mladineo2015). Indeed, the last authors have found that a significant sub-structuring was reported between mtDNA cox2 sequences data sets of the western and eastern populations of A. pegreffii of the Pacific Ocean, obtained from sequences deposited in Genbank, with respect to population of the same species collected by those authors from the Adriatic Sea. On the contrary, the same authors (Blažeković, et al., Reference Blažeković, Pleić, Đuras, Gomerčić and Mladineo2015) did not report significant genetic sub-structuring between Adriatic Sea and other Mediterranean populations of A. pegreffii; this last finding seems to be congruent with the low level of genetic differentiation here estimated at the SSRs level, between the populations of the Adriatic and Tyrrhenian Sea waters. The finding of a genetic substructuring of A. pegreffii from Austral and Boreal regions suggests that geographic distance would restrict the gene flow between ‘antipodean' populations of this parasite species. However, based on the F st values observed at the five SSRs loci, the species A. pegreffii showed a high-gene flow value (N m = 7.8) between the analysed Boreal and Austral populations. This resulted to be at similar degree as that previously found at allozymes level (Mattiucci et al., Reference Mattiucci, Nascetti, Cianchi, Paggi, Arduino, Margolis, Brattey, Webb, D'Amelio, Orecchia and Bullini1997; Mattiucci and Nascetti, Reference Mattiucci and Nascetti2008). Different populations and species of cetaceans, among those analysed in the present study, would be responsible for shaping the observed genetic sub-structuring of Boreal and Austral subpopulations of A. pegreffii, as also suggested in our previous studies (Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014), and by other authors (Blažeković et al., Reference Blažeković, Pleić, Đuras, Gomerčić and Mladineo2015). The low Nei's D and F st values observed between the studied populations of A. simplex (s. s.), revealed a low level of genetic sub-structuring among the analysed samples collected from different basin waters of the NE Atlantic Ocean. This finding was in accordance with previous observations obtained from similar conspecific populations of A. simplex (s. s.) analysed at other nuclear markers (allozymes) (Mattiucci et al., Reference Mattiucci, Nascetti, Cianchi, Paggi, Arduino, Margolis, Brattey, Webb, D'Amelio, Orecchia and Bullini1997), and at mitochondrial gene level (mtDNA cox2) (Mattiucci et al., Reference Mattiucci, Giulietti, Paoletti, Cipriani, Gay, Levsen, Paoletti and Nascetti2018b).
Further, the F IT value higher than zero, as found in both A. pegreffii and A. simplex (s. s.) (Suppl. Table 1a and b), seems to indicate that a certain subdivision between infrapopulations of parasites from different hosts would exist. However, the high level of polymorphism observed at the scored SSRs loci, would require a larger number of specimens to be studied in all the parasite populations collected from different host species, in order to find causes of those differences at the infrapopulation level.
Finally, higher estimates of Nei's D genetic distance was inferred at the interspecific level between A. pegreffii and A. simplex (s. s.) (on average, D Nei = 0.91), based on five considered SSRs loci, with respect to that previously observed between the two sibling species, as estimated from allozyme loci (on average, D Nei = 0.36). However, the latter value was based on genetic data inferred from a larger number of loci scored (N = 24) in populations of the two Anisakis species (Mattiucci et al., Reference Mattiucci, Nascetti, Cianchi, Paggi, Arduino, Margolis, Brattey, Webb, D'Amelio, Orecchia and Bullini1997, Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014).
Concluding remarks
Among the nuclear codominant markers applied to discriminate the two species (i.e. A. simplex (s. s.) and A. pegreffii) of the A. simplex (s. l.) complex, so far allozymes and EF1 α-1 nDNA had provided valuable genetic data for the distinctiveness, reproductive isolation, and absence of gene flow between A. simplex (s. s.) and A. pegreffii (Mattiucci et al., Reference Mattiucci, Cipriani, Webb, Paoletti, Marcer, Bellisario, Gibson and Nascetti2014, Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a, Reference Mattiucci, Giulietti, Paoletti, Cipriani, Gay, Levsen, Paoletti and Nascetti2018b). In the present study, a new set of nuclear markers (SSRs markers) was selected and validated in the two Anisakis species. These novel markers clearly permit to further confirm that A. simplex (s. s.) and A. pegreffii are two distinct species. This finding also remarks the utility of the detected SSRs in the identification of the two species, since SSRs achieve a high-discriminatory power in a nuclear multilocus genotyping approach.
Moreover, the latter approach would provide a powerful mean to investigate, in future researches, various microevolutionary aspects in anisakid taxa, as well as the detection of possible patterns of hybridization/introgression events between the two species in sympatric areas (Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a; and pers.com.). Thus, the developed SSRs loci will be applicable for studies of A. simplex (s. s.) and A. pegreffii from other globally distributed oceanographic waters, where they are infecting several other host species. Furthermore, they will also permit to infer data on the genetic structure of other genetically closely related species, such as A. berlandi, through cross-amplification studies.
SSRs markers have enabled to investigate the genetic variability and population genetic structure of A. pegreffii and A. simplex (s. s.). High genetic variability and low level of population genetic structure was generally found in the two species. This finding would be the result of a large effective population size generally observed in these parasites species, at both intermediate/paratenic (fish) and definitive (cetaceans) hosts level, and the high level of gene flow between populations of the two species, maintained by the large migration habits of some hosts through geographically distant basin waters. The differentiation observed among populations of the two species was generally low. Values found among Boreal and Austral populations within the species A. pegreffii could be explained by smaller population size which could characterize the population of A. pegreffii from the Mediterranean Sea, with respect to that from the Austral Region; however, further genetic analysis and data on the parasite population density of this parasite species are required to support this hypothesis.
Finally, the SSRs markers could be potentially used for wider and advanced population genetics analyses, to investigate a genetic background of parasite adaptation to different host species, and to monitor the genetic variability values of populations of those parasite species. This last aspect would represent a further tool to investigate the possible use of Anisakis species and their genetic diversity at nuclear level, as ecological indicators of the impact of anthropogenic events and of global climate changes, providing insights on the stability of trophic webs of marine ecosystems, as previously demonstrated in anisakids on the basis of other nuclear markers (Mattiucci and Nascetti, Reference Mattiucci and Nascetti2007, Reference Mattiucci and Nascetti2008).
As it has been recently outlined (Cole and Viney, Reference Cole and Viney2018), new methods of analysing population genetics of parasitic nematodes will be available in future studies, such as those gathered from whole-genome sequencing analysis and (dd)RADSeq. Thus, our understanding of nematodes population genetics, biology, evolutionary history and how they are able to respond to selective pressure will dramatically increase. This will likely help to provide a wider comprehensive picture of the ecological and epidemiological drivers (Mattiucci et al., Reference Mattiucci, Cipriani, Levsen, Paoletti and Nascetti2018a) shaping the population genetic structure of the species of the A. simplex (s. l.) complex.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S003118201900074X
Author ORCIDs
Simonetta Mattiucci, 0000-0003-3881-095X.
Acknowledgements
The Authors are very grateful to the PhD students, Drs Marialetizia Palomba and Lucilla Giulietti, for their technical support in the collection of the Anisakis spp. samples from fish hosts, and in the parasites' DNA extraction.
Financial support
SM, MP, AL, PC and GN carried out part of the research with funding received from the ‘European Union Seventh Framework Programme for research, technological development and demonstration', under the grant agreement no. 312068 ‘PARASITE’.
Additionally, SM carried out part of the research work with a grant from “Istituto Pasteur Italia-Fondazione Cenci Bolognetti, Italy”.
Conflict of interest
None.
Ethical standards
We declare that there are no ethical issues related to this paper. We have not used laboratory animals.
We confirm that for vertebrate animals (dead animals) examined for parasites, we have followed the standard principles defined by the European Convention for the Protection of Vertebrate Animals.











