Introduction
Avian schistosomes of the genus Austrobilharzia Johnston, 1917 (Trematoda: Schistosomatidae) are cosmopolitan parasites infecting a wide range of shorebirds in Australia, Canada, United States and Hawaii (Barber & Caira, Reference Barber and Caira1995). Their wide range of distribution can be attributed to migratory birds, mainly seagulls that serve as their definitive hosts (Appleton, Reference Appleton1983; Abdul-Salam & Sreelatha, Reference Abdul-Salam and Sreelatha2004). These schistosomes also use several species of marine snails as intermediate hosts, resulting in a very broad geographical distribution. These hosts include Cerithidea californica, C. scalariformis and Ilyanassa obsoleta from North America (Miller & Northup, Reference Miller and Northup1926; Penner, Reference Penner1950; Short & Holliman, Reference Short and Holliman1961; Barber & Caira, Reference Barber and Caira1995; Leighton et al., Reference Leighton, Ratzlaff, McDougall, Stewart, Nadan and Gustafson2004), C. costata from the Caribbean (Cable, Reference Cable1956), Littorina pintado from Hawaii (Chu & Cuttress, Reference Chu and Cuttress1954), Planaxis sulcatus and Pyrazus australis from Australia (Bearup, Reference Bearup1955; Appleton, Reference Appleton1983), and P. sulcatus and C. cingulata from Kuwait (Abdul-Salam & Sreelatha, Reference Abdul-Salam and Sreelatha1998, Reference Abdul-Salam and Sreelatha2004). The adult flukes usually inhabit the bird's mesenteric veins (Chu & Cuttress, Reference Chu and Cuttress1954; Appleton, Reference Appleton1983) but their cercariae draw attention mainly because of their effect on human health (Gentile et al., Reference Gentile, Picot, Bourdeau, Bardet, Kerjan, Piriou, Le Guennic, Bayssade-Dufour, Chabasse and Mott1996), manifesting as aetiological agents of marine cercarial dermatitis (swimmer's itch) as reported in North America (Grodhaus & Keh, Reference Grodhaus and Keh1958; Sindermann, Reference Sindermann1960; Leighton et al., Reference Leighton, Ratzlaff, McDougall, Stewart, Nadan and Gustafson2004), Australia (Bearup, Reference Bearup1955; Rohde, Reference Rohde1977), Hawaii (Arnold & Bonnet, Reference Arnold and Bonnet1950; Chu, Reference Chu1952; Chu & Cuttress, Reference Chu and Cuttress1954) and China (Liu & Pai, Reference Liu and Pai1974).
The life cycle and the surface topography of the larval stages of several Austrobilharzia species have been described. These included A. variglandis (Chu & Cuttress, Reference Chu and Cuttress1954), A. terrigalensis (Bearup, Reference Bearup1955; Rohde, Reference Rohde1977) and A. penneri (Holliman, Reference Holliman1961). In Kuwait, Abdul-Salam & Sreelatha (Reference Abdul-Salam and Sreelatha1998, Reference Abdul-Salam and Sreelatha2004) reported schistosome cercariae from two marine gastropods, C. cingulata (Gastropoda: Potamididae) and P. sulcatus (Gastropoda: Planaxidae), both with a very low prevalence of infection (less than 1%). The schistosome cercariae in Kuwait have been identified to the genus level as Austrobilharzia sp. based on their surface microtopography, while their species remained unresolved (Abdul-Salam & Sreelatha, Reference Abdul-Salam and Sreelatha2004). It is often difficult to distinguish cercariae of different avian schistosomes based on morphological characteristics alone (Hertel et al., Reference Hertel, Hamburger, Haberl and Haas2002; Aldhoun et al., Reference Aldhoun, Kolárová, Horak and Skirnisson2009; Brant & Loker, Reference Brant and Loker2009). Therefore, the application of molecular techniques is necessary for the identification of the avian schistosome causing dermatitis, since cercariae of different avian schistosomes (Trichobilharzia, Austrobilharzia, Gigantobilharzia) were reported as causative agents of the disease (Kolárová, Reference Kolárová2007). Sequence analysis of nuclear ribosomal DNA and mitochondrial genes is widely used in the identification of schistosome cercariae. This is attributed to the relative identity of 28S, 18S, internal transcribed spacer (ITS) regions and mitochondrial cytochrome oxidase I (mtCO1) across different stages of the digenean life cycle and their ability to differentiate reliably between taxa (Kane & Rollinson, Reference Kane and Rollinson1994; Barber et al., Reference Barber, Mkoji and Loker2000; Dvorák et al., Reference Dvorák, Vanacova, Hampl, Flegr and Horák2002; Lockyer et al., Reference Lockyer, Olson, Østergaard, Rollinson, Johnston, Attwood, Southgate, Horak, Snyder, Le, Agatsuma, McManus, Carmichael, Naem and Littlewood2003; Snyder, Reference Snyder2004; Nolan & Cribb, Reference Nolan and Cribb2005; Vilas et al., Reference Vilas, Criscione and Blouin2005; Brant et al., Reference Brant, Morgan, Mkoji, Snyder, Rajapakse and Loker2006; Aldhoun et al., Reference Aldhoun, Kolárová, Horak and Skirnisson2009; Brant & Loker, Reference Brant and Loker2009; Wang et al., Reference Wang, Li, Ni, Zhai, Chen, Chen and Zhu2009; Karamian et al., Reference Karamian, Aldhoun, Maraghi, Hatam, Farhangmehr and Sadjjadi2011).
Since the schistosome cercariae found in Kuwait Bay had been previously identified as Austrobilharzia species (Abdul-Salam & Sreelatha, Reference Abdul-Salam and Sreelatha2004), in the present study, we aimed to identify the schistosome cercariae using molecular techniques. Sequences of 28S, 18S, ITS2 and mtCO1 regions of the samples were used to identify the species and the phylogenetic relationship of the schistosome species from Kuwait Bay with other schistosomes whose sequences have been deposited in GenBank.
Materials and methods
Sample collection
The snail Cerithidea cingulata (Gastropoda: Potamididae) was collected at low tide from Shuwaikh Bay (29°22'N, 47°58'E) covering 2 km west of Kuwait City. In total, 1475 snails were hand picked from mud between April 2009 and December 2009. In the laboratory, snails were examined for the presence of sporocysts using a dissecting microscope. The identification of larval digenea was made using 0.5 neutral red vital stain specimens (sporocysts and cercariae) according to Schell (Reference Schell1985) and Cable (Reference Cable1956). For each schistosomatidal infection, the sporocysts were cleaned from the snail tissues and stored in 95% ethanol for subsequent DNA extraction.
DNA extraction and polymerase chain reaction amplification
About 50 sporocysts per infected snail were washed with milli-Q filtered water (mqH2O) to remove traces of ethanol. Total gDNA was extracted using proteinase K digestion and organic/aqueous phase extraction (Barber et al., Reference Barber, Mkoji and Loker2000; Al-Kandari & Al-Bustan, Reference Al-Kandari and Al-Bustan2010).
Polymerase chain reactions (PCRs) were performed in a reaction volume of 25 μl (2.5 μl of 10 × PCR buffer, 2.5 mm MgCl2, 10 μm of each deoxynucleoside triphosphate (dNTP), 1.5 μm of each primer and 1 Unit of Amplitaq DNA polymerase (Applied Biosystems, Foster City, California, USA)). The 28S partial DNA sequence was amplified in three overlapping sections using a combination of the primers U178, L1642, U1148, L2450, U1846 and L3449 (Lockyer et al., Reference Lockyer, Olson, Østergaard, Rollinson, Johnston, Attwood, Southgate, Horak, Snyder, Le, Agatsuma, McManus, Carmichael, Naem and Littlewood2003). The PCR conditions consisted of an initial denaturation step at 94°C for 5 min followed by 40 cycles of 30 s at 94°C, 30 s at 52°C and 2 min at 72°C, followed by a final extension of 7 min at 72°C and maintained at 4°C.
The 18S partial DNA sequence was amplified in two overlapping sections using the primer combination primer 1+ primer 2 and primer 3+ primer 4 (Barker et al., Reference Barker, Blair, Cribb and Tonion1993). The PCR conditions used consisted of an initial denaturation step at 95°C for 3 min followed by 35 cycles of 30 s at 95°C, 30 s at 64°C and 2 min at 72°C, followed by a final extension of 7 min at 72°C and maintained at 4°C.
The ITS2 region of the nuclear rDNA was amplified by PCR from the parasite gDNA using two different platyhelminth-specific primers (Eurogenetec, Belgium), its3'trem and its4'trem (Dvorák et al., Reference Dvorák, Vanacova, Hampl, Flegr and Horák2002). The amplification profile consisted of an initial denaturation step at 95°C for 10 min followed by 35 cycles of 30 s at 95°C, 30 s at 65°C and 3 min at 72°C with a final extension of 10 min at 72°C and refrigerated at 4°C. For the amplification of mtCO1, two different trematode specific primers were used, JB3' (Bowles et al., Reference Bowles, Blair and McManus1995) and CO1-R trema (Miura et al., Reference Miura, Kuris, Torchin, Hechinger, Dunham and Chiba2005). The amplification profile consisted of an initial denaturation step at 94°C for 1 min followed by 35 cycles of 30 s at 94°C, 15 s at 46.5°C and 1 min at 72°C with a final extension of 7 min at 72°C and refrigerated at 4°C.
PCR products were excised from an agarose gel and purified using a Nucleospin Extract II kit (Macherey Nagel, Germany) according to the manufacturer's protocol.
DNA sequencing
The purified PCR products were then sequenced in the Gene Analyzer 3130XL using the Big Dye Terminator kit (Version 3.1) according to the manufacturer's instructions (Applied Biosystems). The same amplification primers for the 28S, 18S, ITS2 and mtCO1 were used for the forward and reverse sequencing in separate reactions. Sequences were analysed using EMBOSS (Rice et al., Reference Rice, Longden and Bleasby2000) and aligned using ClustalW2 (Thompson et al., Reference Thompson, Higgins and Gibson1994). The sequences obtained in this study for 28S, 18S, ITS2 and mtCO1 were deposited in GenBank under the accession numbers JF742195, JF742194, HQ106460 and HQ106461, respectively.
Molecular phylogenetic analysis
The sequence data were aligned with other trematode sequences from the family Schistosomatidae deposited in GenBank, using ClustalW2 (Thompson et al., Reference Thompson, Higgins and Gibson1994). Phylogenetic analyses using Maximum Likelihood (ML) and Minimum Evolution (ME) were carried out using MEGA 5.0 (Kumar et al., Reference Kumar, Dudley, Nei and Tamura2008) and Bayesian probabilities using MR BAYES (Heulsenbeck & Ronquist, Reference Heulsenbeck and Ronquist2001). Maximum Likelihood calculations were performed using the Tamurai-nei method as nucleotide substitution model with a uniform substitution rate. The Minimum Evolution analysis was performed using a close-neighbour-interchange algorithm with search level 1, where the evolutionary distances were calculated using the Number of Differences method. The neighbour-joining algorithm (Saitou & Nei, Reference Saitou and Nei1987) was used to generate the initial tree. Bayesian inference analyses were run using MR BAYES (Heulsenbeck & Ronquist, Reference Heulsenbeck and Ronquist2001) employing a general-time-reversible (GTR) model for estimation of invariant (I) and gamma distribution (γ) among site rate variation. Posterior probabilities were estimated for 1,000,000 generations, sampling every 100 generations. The first 250 trees were discarded as burn-in. In all phylogenetic analyses Spirorchis scripta and Spirorchid sp. (Digenea: Spirorchidae) were used as the outgroup. The values of 70% and above in the bootstrap test of phylogenetic accuracy indicated reliable grouping among different members of the family Schistosomatidae.
Results
Frequency of infection
Of the 1475 specimens of Cerithidea cingulata examined for larval digenea, 579 (39.3%) were infected with different larval digeneans. Of the total digenean infections, only five snails (0.9%) were infected with the schistosome sporocysts.
PCR amplification and sequence analysis
The expected PCR products were successfully obtained from sporocysts using the primers stated above for 28S, 18S, ITS2 and mtCO1 of sporocysts. The amplified 28S and 18S regions represented partial sequences and were found to be 3813 bp and 2062 bp, respectively. The amplified ITS2 rDNA complete region and mtCO1 partial region were estimated to be 309 bp and 819 bp, respectively. The obtained sequences were compared among all sporocyst samples and were found to be identical in length and composition.
The obtained sequences of 28S, 18S, ITS2 and mtCO1 were subjected to homology sequence similarity search with other trematode sequences in GenBank (table 1). Analysis of the obtained 28S partial sequence by BLAST search revealed 98.6% identity with A. terrigalensis and 98.5% with A. variglandis. Moreover, analysis of the obtained 18S partial sequence by BLAST revealed 99% identity with A. terrigalensis and 98.3% with A. variglandis (table 2).
Table 1 Digenean species including Austrobilharzia sp. from Kuwait Bay used for phylogenetic analysis and their respective GenBank accession numbers; all species belong to the family Schistosomatidae except for two species in the Spirorchidae (*).

Table 2 A comparison between 28S, 18S, ITS2, mtCO1 sequences of Austrobilharzia sp. from Kuwait Bay and various avian schistosome species deposited in GenBank.

In case of the ITS2 complete sequence, comparison with A. terrigalensis and A. variglandis was not possible because of the absence of the ITS2 sequences in GenBank (table 1). BLAST searches showed that the obtained mtCO1 partial sequence, which codes for the haem-copper oxidase subunit I, had 89% identity with A. terrigalensis and 90% identity with A. variglandis (table 2).
Phylogenetic analysis
The phylogenetic analysis was performed by comparing the obtained 28S, 18S, ITS2 and mtCO1 sequences with other sequences for members of the family Schistosomatidae available in GenBank. Phylogenetic analysis of the Austrobilharzia sp. from Kuwait Bay based on the 28S and 18S partial sequences using various distance and discrete methods (Maximum Likelihood, Minimum Evolution and Bayesian Inference) revealed that the Austrobilharzia sp. from Kuwait Bay clusters with A. terrigalensis and A. variglandis, forming a well-supported monophyletic clade (figs 1 and 2). These results were further supported by the phylogenetic trees constructed based upon the combined dataset of 18S-28S-mtCO1 partial sequence which showed a close relationship of the Austrobilharzia from Kuwait Bay with A. terrigalensis and A. variglandis (fig. 3). Bootstrapping of the sequences with Maximum Likelihood, Minimum Evolution and Bayesian Inference (98, 90, 100 respectively) revealed further significant support for the clade containing A. terrigalensis, A. variglandis and Austrobilharzia sp. from Kuwait Bay. The phylogenetic trees were built based on 28S, 18S and mtCO1 sequences only because of the lack of ITS2 sequences of other Austrobilharzia sp. in GenBank.

Fig. 1 Maximum likelihood tree based on the 28S sequences of Austrobilharzia sp. from Kuwait Bay and other digenean species available in GenBank. Node support is indicated by Maximum Likelihood and Minimum Evolution bootstrap values as well as Bayesian posterior probabilities, respectively. The asterisk indicates 100% nodal support and the dash indicates no significant node support. Branch support is shown only for the major clades.

Fig. 2 Maximum likelihood tree based on the 18S sequences of Austrobilharzia sp. from Kuwait Bay and other digenean species available in GenBank. Node support is indicated by Maximum Likelihood and Minimum Evolution bootstrap values as well as Bayesian posterior probabilities, respectively. The asterisk indicates 100% nodal support and the dash indicates no significant node support. Branch support is shown only for the major clades.

Fig. 3 Maximum likelihood tree based on the combined 18S-28S-mtCO1 sequences of Austrobilharzia sp. from Kuwait Bay and other digenean species available in GenBank. Node support is indicated by Maximum Likelihood and Minimum Evolution bootstrap values as well as Bayesian posterior probabilities, respectively. The asterisk indicates 100% nodal support and the dash indicates no significant node support. Branch support is shown only for the major clades.
Discussion
As the results indicate, the low frequency of Austrobilharzia sp. infection in the snail C. cingulata at Kuwait Bay is consistent with the findings on other avian schistosomes (Grodhaus & Keh, Reference Grodhaus and Keh1958; Loy & Haas, Reference Loy and Haas2001; Karamian et al., Reference Karamian, Aldhoun, Maraghi, Hatam, Farhangmehr and Sadjjadi2011). The observed infection frequency in snails at Kuwait Bay could be related to the migratory behaviour of birds that act as the definitive host in the life cycle of the parasite, resulting in seasonal fluctuations in Austrobilharzia infections in snails. More than 40 species of migratory birds coming from the north pass across Kuwait Bay during the short months of spring and autumn (Haynes, Reference Haynes1979). It is also possible that the observed infection frequency is related to the sheer number of snail intermediate host and the low probability of finding snails infected with Austrobilharzia sp. among the large density of the snail population.
In the present study, we characterized Austrobilharzia from Kuwait Bay using 28S, 18S, ITS2 and mtCO1 sequences. The complete sequence of these genetic markers match each other among all sporocyst samples obtained in this study, indicating the presence of a single schistosome species. It is also worth noting that the 28S and 18S sequences showed no intraspecific variations. The high identity in the 28S and 18S sequences (98.5% and 98.3% respectively) among the schistosomes from Kuwait Bay and other Austrobilharzia species in GenBank confirmed the morphologically based generic identification. On the other hand, the divergence in 28S and 18S sequences between Austrobilharzia sp. from Kuwait Bay and A. terrigalensis and A. variglandis greatly exceeded the intraspecific variation in the sequences between trematode species, suggesting a different species. These results are further supported by a phylogenetic analysis where Austrobilharzia sp. from Kuwait Bay clustered with A. terrigalensis and A. variglandis in a well-supported monophyletic clade. Moreover, the different Austrobilharzia species grouped with Ornithobilharzia in a clade, as previously described by Snyder & Loker (Reference Snyder and Loker2000), Lockyer et al. (Reference Lockyer, Olson, Østergaard, Rollinson, Johnston, Attwood, Southgate, Horak, Snyder, Le, Agatsuma, McManus, Carmichael, Naem and Littlewood2003), Brant et al. (Reference Brant, Morgan, Mkoji, Snyder, Rajapakse and Loker2006) and Brant & Loker (Reference Brant and Loker2009).
Initially, the ITS2 region was used to identify Austrobilharzia sp. from Kuwait Bay to the species level. However, sequence comparison of Austrobilharzia sp. from Kuwait Bay with other Austrobilharzia species using ITS2 sequences was not possible due to the absence of such sequences in GenBank. The ITS sequences of other Austrobilharzia sp. should be included in future studies in order to resolve the relationship and phylogeny among different Austrobilharzia species and validate the use of ITS2 in Austrobilharzia species identification.
Since mtCO1 was shown to be an accurate marker in molecular studies of helminth taxonomy (Vilas et al., Reference Vilas, Criscione and Blouin2005), sequence analysis of the mtCO1 was performed in an attempt to reveal the specific identity of Austrobilharzia from Kuwait Bay. The analysis revealed high sequence variation (10–11%) between Austrobilharzia sp. from Kuwait Bay and Austrobilharzia spp. in GenBank. The same pattern of variation was observed among avian schistosomes of the genera Austrobilharzia and Trichobilharzia (Brant & Loker, Reference Brant and Loker2009) where species exhibited mtCO1 sequence variation of 13.4%. Therefore, Austrobilharzia sp. from Kuwait Bay belongs to neither A. terrigalensis nor to A. variglandis, and is likely to be an Austrobilharzia species that has not been characterized at the molecular level. Despite the high interspecific differences in the mtCO1 sequences among Austrobilharzia species, Austrobilharzia sp. from Kuwait Bay clustered with A. terrigalensis and A. variglandis forming a monophyletic clade in the phylogenetic tree constructed based on combined dataset of 18S-28S-mtCO1 sequences. However, the cercariae of Austrobilharzia sp. from Kuwait Bay may be more closely related to other Austrobilharzia species, e.g. A. penneri, that has been characterized morphologically (Short & Holliman, Reference Short and Holliman1961) but not molecularly. The Austrobilharzia from Kuwait Bay and A. penneri not only are morphologically similar but both also successfully infect snails of the genus Cerithidea (Short & Holliman, Reference Short and Holliman1961). The molecular relationship between the two species could not be investigated further due to the lack of DNA sequences of A. penneri in GenBank.
Results from the present study support previous molecular studies performed on the family Schistosomatidae. The topology of the Schistosomatidae reported by the 28S and the combined dataset of 18S-28S-mtCO1 sequences was similar to that produced by Snyder (Reference Snyder2004) and Brant et al. (Reference Brant, Morgan, Mkoji, Snyder, Rajapakse and Loker2006). All three studies reported the existence of an Ornithobilharzia + Austrobilharzia clade as basal to the rest of the schistosomes. This was not supported by the phylogenetic tree based on the 18S sequences, where the Heterobilharzia + Schistosomatium clade falls basal to the rest of the schistosomes, a topology that matches that of Snyder (Reference Snyder2004). The difference in the 18S tree topology could be reflecting the evolution of the gene and not necessarily the evolutionary history of the species (Littlewood & Olson, Reference Littlewood and Olson2001).
The present study confirmed the morphological identification of Kuwait's schistosome cercariae as Austrobilharzia sp. Its accurate species identification and classification requires that molecular identification of Austrobilharzia species from other continents be expanded. Only then can the accurate identification and relationship of different species become clearer.
Acknowledgements
The authors extend their deepest appreciation and gratitude to Dr Jitka Aldhoun of the Department of Zoology, The Natural History Museum, London, for assistance with the PCR primers. Authors also gratefully acknowledge Professor Armand Kuris of the Department of Ecology, Evolution and Marine Biology, University of California, Santa Barbara, USA for reviewing the manuscript. Molecular analysis was performed with the support of the General Facility Project (GS 01/02) in the use of the ABI 3130xl Gene Analyzer at the Biotechnology Center, Faculty of Science, Kuwait University.