Loads of trematodes: discovering hidden diversity of paramphistomoids in Kenyan ruminants

Summary

Paramphistomoids are ubiquitous and widespread digeneans that infect a diverse range of definitive hosts, being particularly speciose in ruminants. We collected adult worms from cattle, goats and sheep from slaughterhouses, and cercariae from freshwater snails from ten localities in Central and West Kenya. We sequenced cox1 (690 bp) and internal transcribed region 2 (ITS2) (385 bp) genes from a small piece of 79 different adult worms and stained and mounted the remaining worm bodies for comparisons with available descriptions. We also sequenced cox1 and ITS2 from 41 cercariae/rediae samples collected from four different genera of planorbid snails. Combining morphological observations, host use information, genetic distance values and phylogenetic methods, we delineated 16 distinct clades of paramphistomoids. For four of the 16 clades, sequences from adult worms and cercariae/rediae matched, providing an independent assessment for their life cycles. Much work is yet to be done to resolve fully the relationships among paramphistomoids, but some correspondence between sequence- and anatomically based classifications were noted. Paramphistomoids of domestic ruminants provide one of the most abundant sources of parasitic flatworm biomass, and because of the predilection of several species use Bulinus and Biomphalaria snail hosts, have interesting linkages with the biology of animal and human schistosomes to in Africa.

INTRODUCTION

The Superfamily Paramphistomoidea is a prominent group of digeneans where adults are characterized by the absence of an oral sucker and the presence of an acetabulum at or near the posterior end of the body. The systematics of this group of digeneans is a work in progress. Sey (1991) concluded it is comprised of eight families, whereas Jones (2005a ) concluded there are 12 families. Paramphistomoids are often called rumen flukes because many of the best-known representatives live in this habitat in domestic ruminants. However, many species also inhabit the intestines of fish, amphibians, reptiles, birds and non-ruminant mammals. They feature a life cycle in which cercariae produced in rediae emerge from snails and encyst on vegetation as metacercariae, which are later ingested by the definitive host (Jones, 2005a ). As part of a larger study to determine how digenean community diversity influences the transmission of schistosomes in Kenya, we provide new results regarding the overall diversity and host relationships of paramphistomoids in Kenya, based on cercariae collected from snails and adult worms from domestic animals from abattoirs.

Paramphistomoids are of interest to parasitologists in several contexts. They are diverse in number of species and provide an understudied model group for those focused on revealing patterns and mechanisms of diversity. Of the 12 recognized paramphistomoid families recognized by Jones (2005a ), representatives of nine occur in Africa. The diversity of paramphistomoids in Africa reflects the presence of many species of terrestrial mammals, including elephants, rhinoceroses, hippopotami and a rich diversity of wild and domestic ruminants. Three families in particular (Paramphistomidae, Gastrodiscidae and Gastrothylacidae) are speciose in Africa. The distribution of diversity in rumen hosts can partly be explained by characters (e.g. regressed pharyngeal appendages) that are apomorphic, which have allowed them to colonize the forestomach (Sey 1991). The three families comprise over 40% of all known paramphistomoids, the majority of which use ruminants as their definitive hosts (Sey, 1991).

Paramphistomoids have thick bodies, which make detailed morphological characterization of adult features and species identification challenging (Horak, 1971; Jones, 1991; Mage et al. 2002; Rinaldi et al. 2005). The bodies of paramphistomoid cercariae are also relatively thick and typically filled with cystogenous material or pigment, also rendering identification difficult. Nonetheless, a meticulous framework for paramphistomoid identification and classification has been developed (see reviews by Sey, 1991; Jones, 2005a ). Given the inherent difficulties in identification, coupled with a growing list of studies from other digenean groups documenting the presence of cryptic species (Detwiler et al. 2012; Herrmann et al. 2014; McNamara et al. 2014), paramphistomoids are ideal for studies attempting to meld traditional morphological identification with sequence data characterization provided by molecular approaches. The number of studies that use molecular techniques to provide assessments of the diversity of paramphistomoids have in general been limited, especially so for African species (Lotfy et al. 2010; Mansour et al. 2014; Sibula et al. 2014; Titi et al. 2014; Dube et al. 2015).

In addition to being speciose, paramphistomoids are often remarkably abundant (Horak, 1971; Cheruiyot and Wamae, 1988; Rolfe et al. 1994; Sanabria and Romero, 2008). In fact, one might be hard pressed to find a larger source of sheer digenean biomass than is presented routinely at abattoirs by ruminant paramphistomoids. Given the large worm populations that can occur in individual cattle, goats or sheep, vast numbers of paramphistomoid eggs are regularly passed into the environment. In rural West Kenya, we can routinely collect 10 000 paramphistomoid eggs from a single cow dung sample. As domestic ruminants regularly seek water from natural habitats, it is not surprising that many paramphistomoid eggs enter freshwater, creating the potential for high levels of infection in their snail hosts (Chingwena et al. 2002a ; Mohammed et al. 2016).

A review of the East African paramphistomoid literature reveals that many of the described species are transmitted by Biomphalaria and Bulinus, the snail genera also of concern with respect to their role in transmission of human schistosomiasis in Africa (Dinnik, 1954; Dinnik and Dinnik, 1957; Dinnik, 1961; Eduardo, 1983; Brown, 1994; Chingwena et al. 2002b ; Jones, 2005b , c ). In some areas, Bulinus and Biomphalaria are the most commonly implicated snail hosts for paramphistomoids (Dinnik, 1965; Wright et al. 1979; Loker et al. 1981; Chingwena et al. 2002b ; Ahmed et al. 2006; Mohammed et al. 2016). The presence of other digenean species utilizing the same snail species as schistosomes could be a factor that influences the overall success of animal and human schistosome transmission (Lim and Heyneman, 1972; Combes, 1982; Hechinger et al. 2011; Spatz et al. 2012). This is particularly so for species such as paramphistomoids that produce rediae as larval stages within their snail hosts, because rediae may attack, damage and consume schistosome sporocysts (Lim and Heyneman, 1972).

We collected cercariae and adult worms from ten localities in Kenya. We provide stained whole mounts and provisional identification of adults that are linked to sequence data for cytochrome oxidase 1 (cox1) and the internal transcribed region 2 (ITS2). In some cases, we provide matches with sequences obtained from cercariae and adult worms thus providing probable life cycle linkages. We also provide new hypotheses for phylogenetic relationships among the paramphistomoids that include available sequences from NCBI GenBank, which show that some species of paramphistomoids are geographically widespread throughout Africa. Data presented here will contribute to an increased understanding of the superfamily Paramphistomoidea, including providing greater clarification for how these worms are distributed among hosts, their potential roles if any in causing disease in domestic or wild animals, and their interactions with other digeneans, including schistosomes.

MATERIALS AND METHODS

Sampling

We collected larval and adult paramphistomoids from ten different localities in central and especially western Kenya between 2005 and 2015 (Table 1). All species of field-collected aquatic snails were brought to the laboratory at Kisian, near Kisumu, Kenya. The snails were cleaned and then placed individually into 12-well tissue culture plates in 3 mL of aged tap water. The tissue culture plates were placed in natural light for 2 h to induce shedding of cercariae. Snails shedding cercariae were identified using keys and information in Brown and Kristensen (1989) and Brown (1994), and cercariae were preliminarily identified using keys (Frandsen and Christensen, 1984; Schell, 1985) and by reference to regional monographs (e.g. Fain, 1953). All cercariae designated as paramphistomoids were confirmed as such according to Sey (1991). Snails were either dissected at the time of collection to procure rediae, or re-shed two and four weeks later to determine if snails were harboring prepatent infections at the time of collection. Snails were kept in 20 L plastic tanks and fed red leaf lettuce following collection. Cercariae and rediae were preserved in 95% ethanol for later molecular analysis.

Table 1. Collection localities in central and west Kenya

Adults were collected from the rumen or reticulum of Bos indicus, Capra aegagrus hircus and Ovis aries from one slaughterhouse in central Kenya and three in Western Kenya (Table 1). Adults were preserved in 95% ethanol for later molecular and morphological identification.

Staining adult worms

Adult worms were placed into 70% ethanol for 24 h prior to staining. Sections of the adult worms were stained and mounted according to Eduardo (1982). Because of their thickness, each adult was sectioned frontally using a razor blade. Part of the posteroterminally placed acetabulum was severed and used for molecular analysis.

Collection of molecular data

A partial sequence of cox1 mtDNA and internal transcribed spacer two (ITS2) were amplified by polymerase chain reaction (PCR) to facilitate differentiation among paramphistomoid specimens. One to six cercariae, one to three rediae or a portion of the acetabulum from adults were used for DNA extraction. Genomic DNA was extracted from 120 paramphistomoid samples (Table 2) by the alkaline-lysis (HOT-SHOT) method (Truett et al. 2000), or by the QIAamp DNA Micro Kit following the manufacturer's instructions, with a final elution volume of 30 µL (Qiagen, Valencia, CA). Although not the equal of the QIAamp Kit with respect to absolute quality of the DNA produced, the HOT-SHOT method also produced DNA of quality and proved more amenable for use under conditions where controlled conditions were less available.

Table 2. Specimen name, host collected from, collection locality, provisional identification, Museum of Southwestern Biology/KEMRI voucher numbers, and GenBank accession numbers of paramphistomoid specimens used in this study

PA1-PA44 contain representatives of the 16 different clades used to construct the ML and Bayesian trees. PA45-PA120 were included in the preliminary trees. An (*) denotes samples that are in Kenya.

Cox1 oligonucleotide primers were designed based on the barcode region (Folmer et al. 1994) and on conserved regions in the Fasciola hepatica (NC_002546), Paragonimus westermani (AF219379) and Paramphistomum cervi (NC_023095) mitochondrial genomes. Cox1 was amplified using primers 123F [5′-ATTCGTTTGAACTATATGGA-3′] and 858R [5′-CATATGATGAGCCCAAACAAC-3′]. The volume of each PCR reaction was 25 µL with 1 µL of 100 ng of DNA, 0·8 mm L⁻¹ dNTPs, 2·5 mm L⁻¹ MgCl2, 0·25 units of Ex Taq DNA (Clontech, Mountain View, CA) and 0·4 µ m L of each primer. PCR cycles were programmed as follows: 2 min denaturation hold at 94 °C; 94 °C for 1 min, 46 °C for 30 s and 72 °C for 1 min for three cycles; 94 °C for 1 min, 45 °C for 30 s, and 72 °C for 1 min for three cycles; 94 °C for 1 min, 44 °C for 30 s and 72 °C for 1 min for three cycles; 94 °C for 1 min, 44 °C for 30 s and 72 °C for 1 min for 20 cycles, and followed by an extension step for 7 min at 72 °C.

ITS2 was amplified using GA1 [5′-AGA ACA TCG ACA TCT TGA AC-3′] (Anderson and Barker, 1998) and BD2 primers [5′-TAT GCT TAA ATT CAG CGG GT-3′] (Bowles et al. 1995). The volume of each reaction was 25 µL, with 12·5 µL of Premix Taq™ (Clontech, Mountain View, CA), 0·4 µ m L⁻¹ of each primer, and one μL of 55 ng of DNA. PCR cycles were performed on Eppendorf Mastercycler epigradient machines, which were programmed as follows: 1 C s⁻¹ rate of change, one cycle at 98 °C for 10 s, followed by 30 cycles of 98 °C for 1 min, 52 °C for 2 min and 72 °C for 1 min with an extension step for 7 min at 72 °C.

PCR fragments were separated by agarose gel electrophoresis and visualized with 0·5% GelRed™ Nucleic acid gel stain (Biotium, Hayward, CA). PCR products were purified using the QIAquick purification kit (Qiagen, Valencia, CA) or by ExoSap-IT^® (Affymetrix, Santa Clara, CA). Both strands were sequenced using an Applied Biosystems 3130 automated sequencer and BigDye terminator cycle sequencing kit Version 3.1 (Applied Biosystems, Foster City, CA). DNA sequences were verified by aligning reads from the 5′ and 3′ directions using Sequencher 5·0 and manually corrected for ambiguous base calls (Gene Codes, Ann Arbor, MI).

Outgroup determination

To determine the most appropriate outgroup available for our data, we reconstructed trees with the most likely outgroups based on Lockyer et al. (2003) and chose the sister group to the paramphistomoids (ingroup). Species from the following nine families were used from 12 digenean mitochondrial genomes for maximum-likelihood (ML) analysis: Dicrocoelium dendriticum (NC_025280), Fasciola gigantica (NC_024025), P. cervi (NC_023095), Opisthorchis felineus (NC_011127), Clonorchis sinensis (NC_012147), Orthocoelium streptocoelium (NC_028071), Echinostoma hortense (NC_028010), Fischoederius elgonatus (NC_028001), P. westermani (NC_027673), Eurytrema pancreaticum (NC_026916), F. hepatica (NC_002546) and Ogmocotyle sikae (NC_027112).

Sequence alignment and phylogenetic analyses

Phylogenetic analyses were done with cox1 and ITS2 sequences using ML and Bayesian interference (BI). The analysis included four specimens from NCBI-GenBank for cox1 and 43 for ITS2 (Table 2). Non-identical haplotypes of cox1 and ITS2 sequences were aligned by eye and edited in MEGA6 (Tamura et al. 2013). A total of 690 bases were used for cox1 alignment and 385 bases for ITS2 alignments. Sequences generated in this study were submitted to GenBank (Table 2). ML analyses used PAUP* 4·0 b10 (Wilgenbusch and Swofford, 2003) and BI analyses were carried out using MrBayes (v 3.12) (Ronquist and Huelsenbeck, 2003). MrModeltest 2·0 (Nylander, 2004) was used to find the best fit model of substitution for BI and ML for both genes. Heuristic searchers were utilized for ML analyses (excluding the third codon for cox1) and 100 bootstrap replicates were run for each dataset. For BI analyses of the cox1 dataset (excluding the third codon for cox1), the parameters were: nst = 6, rates = invgamma and ngammacat = 4. Four heated chains were run simultaneously for 1 000 000 generations. For BI analyses of the ITS2 dataset, the parameters were: nst = 6, rates = gamma and ngammacat = 4. Four heated chains were run simultaneously for 1 400 000 generations. In both datasets, the trees were sampled every 100 cycles, and the first 25% of trees with pre-asymptotic likelihood scores were discarded as burn-in. A number of generations were determined sufficient because the s.d. dropped below 0·01 at the end of the runs.

Nucleotide substitution saturation at the third codon was tested in DAMBE5 (Xia, 2013) for cox1. Uncorrected pairwise distance values were calculated in MEGA6 (Tamura et al. 2013). Data were summarized within and between groups (Tables 3 and 4). We used similar criteria of other studies that used a P-distance value >5% difference with cox1 and nd1 mtDNA markers and >1·0% for ITS to indicate separate species (Vilas et al. 2005; Brant and Loker, 2009; Detwiler et al. 2010).

Table 3. Intra- and interclade P- distance values of cox1 amplified from paramphistomoids from Kenya

Values in bold are intraclade divergences. Note that “–” indicates only a single specimen was collected and within distances could not be calculated.

Table 4. Intra- and interclade p-distance values of ITS2 amplified from paramphistomoids from Kenya

Values in bold are intraclade divergences. Note that “–” indicates only a single specimen was collected and within distances could not be calculated.

RESULTS

Samples

Paramphistomoid adults were collected from three species of ruminants and cercariae and/or rediae were collected from four different genera of planorbid snails (Biomphalaria, Bulinus, Ceratophallus, Segmentorbis) from ten localities in central and west Kenya (Tables 1 and 2). Paramphistomoid cercariae were not found in other snail species examined (Melanoides tuberculata, Radix natalensis, Physa acuta and Bellamya unicolor). Ruminants were typically heavily infected, and often hundreds of adult worms could be quickly collected per host. From our samples collected, we examined and sequenced 79 adult and 41 cercariae specimens (120 total specimens) that represented obvious variants. To facilitate sampling if a large numbers of adult worms were acquired from a single host, we separated them by differences in adult host morphology (size and presence of a pouch or a genital sucker). To further assure collection of a diversity of specimens, we sampled both adult worms and rediae/cercariae from different localities

Outgroup determination

With the diversity of sequence data available in GenBank, our analysis revealed that O. sikae (Notocotylidae) is more closely related to paramphistomoids than members of Echinostomatidae or Fasciolidae used as outgroups for other paramphistomoid molecular phylogenies (Lotfy et al. 2010; Shylla et al. 2011; Ghatani et al. 2012). For phylogenetic analyses of both genes, we used three species of notocotylids as outgroup taxa.

Cox1 phylogenetic analyses and pairwise distance divergences

In general, trees were first constructed incorporating all 120 specimens (Supplementary Figs. S1 and S2). Because some clades were represented by multiple specimens (haplotypes with a 1–4 bp difference for cox1) we reduced the number of specimens per clade to simplify the trees for display purposes (Figs. 1 and 2). Many of the deeper nodes were not supported; however, the trees nonetheless provided a useful way to visualize the overall diversity of specimens found, and to provide comparisons with available systematic treatments. The specific clades identified (names next to the bolded black vertical lines) on the cox1 tree represent conspecifics (Fig. 1).

Fig. 1. Phylogenetic relationships of 44 samples of paramphistomoids from this study and from GenBank based on cox1 (690 bp) sequences inferred from ML (bootstrap values) analysis. Specimens are named based on sample name, the host it was collected from and are colour coded based on intraclade P-distance values <1·3% and interclade values >6·5%. An orange star represents clades where we matched cercariae and adult sequences. Identifications were made based on GenBank sequences and on the species descriptions in the literature (parentheses). An (*) denotes intermediate host use from studies in the literature that have not been sequenced confirmed.

Fig. 2. Phylogenetic relationships of 44 samples of paramphistomoids from this study and from GenBank based on ITS2 (385 bp) sequences inferred from ML (bootstrap values) analysis. Specimens are named based on sample name, the host it was collected from, and colour coded based on clade designation from cox1 distance values. A red star represents clades where we have found evidence of putative hybrids. Adjacent to these indicated clades, are clade numbers that correspond to the same specimens and clade numbers as appearing on the cox1 tree (Fig. 1).

Partial sequences of cox1 (690 bp) were obtained for all 120 samples (Supplementary Fig. S1). ML and BI (Supplementary Fig. S3) trees were created for the cox1 alignment, and the ML tree is shown (Fig. 1). MrModeltest 2·3 selected the GTR + I + G model of nucleotide substitution. Based on bootstrap and posterior probabilities in Table 3, 16 distinct cox1 clades were identified among Kenyan specimens and are portrayed alongside the tree in Fig. 1 (vertical black lines or arrows). We used genetic distance data to determine if a clade was comprised more than one species. A single species was of determined for specimens with genetic distance values <1·3%, and species were designated as distinct when genetic distance values were >6·2% (Table 3). Most interclade pairwise distance values were >10·0% and they ranged up to 19·9%. These same clade numbers or scientific names were also used adjacent to the ITS2 tree in Fig. 2.

ITS2 phylogenetic analyses and pairwise distance divergences

For ITS2, sequences were obtained from all 120 samples and our phylogenetic analyses also included 46 samples from GenBank (Supplementary Fig. S2). The ITS2 alignment included 61 bp of 5·8S, 283 bp of ITS2 and 46 bp of 28S. The average intraclade pairwise distance was 0·30% and the average interclade pairwise distance was 3·9% (Table 4). MrModeltest 2·3 selected the GTR + G model of nucleotide substitution for ITS2. Both BI and ML analyses were run using 33 or 46, respectively, additional relevant species sequences from GenBank, with the ML tree shown (Figs 2 and Fig S4). Not surprisingly, the degree of resolution provided by phylogenetic analysis of ITS2 sequences was not high given the more conservative rate of change of this widely used nuclear gene marker (Locke et al. 2010). Based on ML and BI analyses, 12 ITS2 clades were identified among our Kenyan specimens (Fig. 2 and Supplementary Fig. S4). Intraclade genetic distance values were <0·6%, and interclade genetic distance values were >1·0%.

Further comparisons of the cox1 and ITS2 datasets

Cox1 and ITS2 trees did not conflict, but the ITS2 trees did not have as much support for the deeper nodes as cox1 (Figs 1 and 2). All 12 clades from ITS2 were represented in the cox1 dataset. The cox1 genetic distance data enabled differentiation among some of the worms clustered with Cotylophoron cotylophorum in the ITS2 dataset, and also clearly differentiated clades 14 and 15 (Fig. 2).

In three cases (clades 4, 10 and 16), cox1 sequence matches (<1·3%) were obtained between worms from ruminants and cercariae from snails (Fig. 1, orange stars). Clade 2 matched an ITS2 sequence from GenBank of cercariae from Ceratophallus natalensis, thus also confirming the intermediate host for this clade (Fig. 1). In four cases (clades 1, 5, 10 and 12), sequences were found from cercariae with no matches from adult worms for either sequence (Fig. 1). In at least five cases (PA7, PA26, PA27, PA35 and PA42), the ITS2 nuclear sequences obtained clustered in different clades than what is seen in the cox1 trees (clades highlighted with red star in Fig. 2). These samples appear to have nuclear mitochondrial discordance (NMD) and are identified as worms with likely hybrid ancestry (see discussion).

Provisional identification of the paramphistomoids

Provisional identifications were based on the paramphistomoid systematics literature (Eduardo, 1983; Sey, 1991; Jones, 2005b , c , d ) pertaining to intermediate or definitive host use, and descriptions of adult worms in comparison to our mounted adult specimens (Table 5, Fig. 3). Some of the sequences we obtained matched sequences from named species in GenBank, and in those cases the names we provide here are the ones from GenBank (clades 4, 8 and 16). Four clades were represented only by cercariae and did not match any sequences derived from adult worms in this study or from GenBank. These included two clades from B. pfeifferi (clades 1 and 12), one from Segmentorbis kanisaensis (clade 5) and one from C. natalensis (clade 10). Our 16 clades represented three different families of Paramphistomoidea: Gastrothylacidae, Paramphistomidae and Stephanopharyngidae. Species names in quotation marks in Fig. 1 were assigned based on our morphological identification from species descriptions.

Fig. 3. Sections of adult paramphistomoids collected from domestic ruminants in Kenya and their provisional identifications. (A) Calicophoron phillerouxi, (B) Calicophoron raja, (C) Calicophoron clavula, (D) Calicophoron microbothrium, (E) Cotylophoron sp., (F) Cotylophoron cotylophorum, (G) Cotylophoron sp., (H) Carmyerius exporous, (I) Carmyerius gregarius, (J) Carmyerius mancupatus. Note that the photographed specimens represent sections of adults, and presence of some organs like the testes (T) or genital sucker (GS) are indicated. For the genus Carmyerius, a ventral pouch was present, but is not visible in the sections chosen for presentation.

Table 5. Provisional identification of the paramphistomoids was based on species descriptions and intermediate host use from the literature and on position in phylogenetic trees

Cercariae (C), adults (A) and their associated hosts are listed. Ventral pouch, acetabulum type and genital sucker were useful morphological features for genus and species placement.

DISCUSSION

Paramphistomoid flukes are speciose in sub-Saharan Africa, reflective of the presence there of many mammal species, particularly wild and domestic ruminants. These flukes are also ubiquitous and can have a high prevalence among domestic ruminants reaching 100% in some villages (Chingwena et al. 2002a ; Nzalawahe et al. 2015). During our sampling of Kenyan slaughterhouses we found up to 90% of the domestic ruminants infected, and many individual animals harboured hundreds of adult worms. Of the many adult worm and cercariae samples collected, we further investigated 120 samples (79 adult worms and 41 cercariae) determined most likely to be genetically distinctive. We found 16 distinct clades in three families of the Paramphistomoidea. For future comparisons, all of our specimens are available as vouchers at the Parasite Division, Museum of Southwestern Biology (MSB) or at the Kenyan Medical Research Institute (KEMRI).

Previous studies have used the easily obtained ITS2 sequence as a molecular marker to distinguish among paramphistomoid species (Itagaki et al. 2003; Rinaldi et al. 2005; Goswami et al. 2009; Lotfy et al. 2010; Sanabria et al. 2011; Ichikawa et al. 2013; Shylla et al. 2013; Ghatani et al. 2014; Dube et al. 2015). ITS2 is helpful for distinguishing paramphistomoid genera and differentiating more divergent species within a genus (Rinaldi et al. 2005; Ghatani et al. 2012). Because mitochondrial DNA accumulates substitutions more frequently than the internal transcribed spacers, it is more useful to differentiate among closely related species, particularly cryptic species (Blouin, 2002; Vilas et al. 2005; Locke et al. 2015), or to reveal intraspecific variation (Ghatani et al. 2014). Consequently, we used genetic distance values for cox1 sequence data as the primary means to delineate species. For cox1, interclade P-distance values were >6·2%, although the majority of pairwise comparisons were >10·0%. In contrast, intraclade pairwise divergence values were <1·3%. Other studies have used a P-distance value >5% difference with cox1 and nd1 mtDNA markers to indicate separate species (Vilas et al. 2005; Brant and Loker, 2009; Detwiler et al. 2010). Our data suggests that ITS2 should not be used alone to differentiate species for paramphistomoids.

We also examined the delineated clades with respect to where they grouped in either ML or BI phylogenetic analyses based on either cox1 or ITS2 sequences. In general, there was low bootstrap/posterior probability support for many of the deeper nodes in either ML or BI trees, suggesting that broader taxon sampling, along with sequencing of additional markers, is needed to more definitively support or refute the morphologically based systematic framework developed for paramphistomoids (Sey, 1991; Jones, 2005a ). The phylogenetic trees were useful, however, in providing preliminary hypotheses for how the various clades were related to one another (see the paragraph below). Relative to other paramphistomoid molecular phylogenetic studies involving specimens from African ruminants and snails, we recovered five out of the six previously reported taxa from Kenya, Egypt and Tanzania noted by Lotfy et al. (2010), three of the three identified taxa from Zimbabwe, Zambia and Botswana (Dube et al. 2015) and one of the two identified taxa from Algeria (Titi et al. 2014). The extent of overlap among specimens recovered from all four studies suggests that at least some of the species have broad distributions in Africa. Additional sampling is needed to provide a more comprehensive picture of African paramphistomoid diversity, particularly from Central and West Africa.

The phylogenetic trees provided support for anatomically based taxon delineations as four clades identified as Calicophoron grouped together, as did three clades of Carmyerius and four clades of Cotylophoron. Furthermore, worms in the Stephanopharyngidae (Stephanopharynx) formed a clade, as did presumptive members of the Gastrothylacidae. However, all presumptive members of the Paramphistomidae did not group together. It is possible that this is a paraphyletic group or certain genera, such as Cotylophoron belong in a different family. Clade 1 is quite divergent from the other specimens discussed and it is possible it represents a different family or superfamily. The trees also show some incongruences between nuclear and mitochondrial sequences (discussed further below).

With respect to host use, specimens from a particular clade were reported from the same snail host species or genus. Also, different clades that group together tend to share the same genus of snail host (Calicophoron, in clades 13–16, in Bulinus) or snail genera in related tribes (Carymerius in clades 2, 3 and 5 in Segmentorbis and Ceratophallus). For 10 of 11 clades for which snail host usage could be identified, those snails belong in the family Planorbidae. Snail host use may thus have had an important impact on paramphistomoid diversification, which has also been suggested for other digenean groups (Brant and Loker, 2013). In only one instance have we found cercariae that we have assigned to the same clade (clade 10) that derive from two different snail genera: cercariae from C. natalensis collected from this study and cercariae from Biomphalaria sudanica collected by Lotfy et al. (2010). Many other digenean groups also indicate high first intermediate host specificity (Shoop, 1988; Donald et al. 2004; Detwiler et al. 2010; Brant and Loker, 2013). By contrast, adult worms of a particular clade were often recovered from more than one definitive host species, and we recovered up to three different taxa of paramphistomoids from an individual bovine.

Sequence data derived from life cycle stages from different hosts provide an important alternative way to piece together the complex life cycles of digeneans, especially when experimental exposures are not possible (Chibwana et al. 2015). We provide supportive evidence for the life cycles of four of our identified clades (Fig. 1) by matching genetic sequences (<0·6% for ITS2 and <1·3% cox1) collected from cercariae and adults: (1) ITS2 sequences from cercariae from C. natalensis (GU735645) collected in Kenya grouped with sequences from adult worms we recovered from cattle (clade 2), provisionally identified as Carmyerius exporous (Dinnik and Dinnik, 1960). (2) Cercariae (clade 4) we collected from C. natalensis matched adults collected in this study as well as two adults from Botswana (KP639636) and Kenya (GU735658) identified as Carmyerius dollfusi by Dube et al. (2015). The latter species was synonymized with C. mancupatus (Sey, 1991), a species known to be transmitted by C. natalensis (Dinnik, 1965). (3) Sequences from seven adults we obtained (clade 15) matched sequences collected from a cercariae sample from B. forskalii. We provisionally identified the adults as C. phillerouxi, which is known to be transmitted by B. forskalii (Dinnik, 1961). (4) Lastly, two cercariae samples we collected from B. forskalii matched with 23 adults collected in this study, and with one cercariae sample from B. forskalii and 18 adults in GenBank, all of which were identified as C. microbothrium (clade 16). As the host record and sequence databases grow, the probabilities that more matches will be found also increases, providing a way forward in working out life cycles that will help offset increasing difficulties in doing so with more conventional experimental infections.

The most common paramphistomoid genus we collected was Calicophoron (40 out of the 120 specimens examined), and the most abundant species was Calicophoron microbothrium which is transmitted by bulinid snails. This species is the most geographically widespread paramphistome in Africa, its presence confirmed with molecular markers from Egypt, Kenya, Tanzania, Zambia, Zimbabwe, South Africa, Algeria and Botswana (Lotfy et al. 2010; Titi et al. 2014; Dube et al. 2015). Given the difficulties in discriminating this species from others based on morphology alone, the broad geographic distribution, and the diversity of different bulinid snails reported as hosts, this species is a good candidate for further inspection as a possible complex of cryptic species. Presently the best sequence available to evaluate this possibility is cox1, but most of the data in the literature thus far for this species are for ITS2. Our ML analysis based on 354 bp of ITS2 (figure not shown) suggests there are distinct clades among the samples identified as C. microbothrium in GenBank, with an average distance among them of 0·75%. Other sequence markers are needed to determine if C. microbothrium is a complex of cryptic species, and how well differentiated they prove to be from the other Calicophoron clades (13–15) identified in this study.

We found some specimens with discordant nuclear and mitochondrial sequences, consistent with the possibility of hybrid origins (red stars, Fig. 2). For example, two samples (PA12 and PA24) grouped with C. microbothrium in the ITS2 trees, but fell in their own clade (3) in the cox1 trees. PA12 and PA24 were also morphologically distinct from C. microbothrium, being provisionally identified as members of the gastrothylacid genus Carmyerius. As we have noted, multiple species of paramphistomoids are frequently recovered from a single ruminant host, creating circumstances conducive for potential hybridization. The putative parental species and hybrids (PA7, PA12, PA24 PA27, PA35) all use Bulinus as intermediate hosts. It seems possible that the likelihood of successful hybridization would be increased if both parental species use the same genus or species of intermediate host, if as appears intermediate host use is more specific than definitive host use among the paramphistomoids. Other examples of sequence discordance in digeneans also involve groups with closely related species that can hybridize, and that share snail hosts, such as with some species of fasciolids and schistosomes (Steinauer et al. 2008; Peng et al. 2009). Further studies using microsatellite markers or RADSeq technology will be needed to verify a hybrid origin for paramphistomoids with discordant sequences.

Members of the basommatophoran family Planorbidae are the most common intermediate hosts transmitting paramphistomoids in Kenya, although snails of the Family Lymnaeidae have also been identified as hosts for paramphistomoids in East Africa (Sey, 1991). The snail hosts for some of the clades we have identified such as clades 3, 6, 7, 8 (C. cotylophorum), 9 and 11 (Stephanopharynx sp.) are unknown or require additional sequence-based verification. Bulinus snails, with an ancient history and diversification in Africa (Van Damme 1984; Brown, 1994; De Groeve, 2005), are particularly prominent as African paramphistomoid hosts (Sey, 1991). By contrast, Biomphalaria supports fewer paramphistomoid species and has a much shorter evolutionary history in Africa, with estimates ranging from <1–5 mya (million years ago) (Woodruff and Mulvey, 1997; Campbell et al. 2000; DeJong et al. 2001). It is noteworthy that clade 1, which is known only from cercariae from B. pfeifferi, is one of the most divergent clades we recovered. Clade 1 cercariae are also much larger than the other paramphistomoid cercariae we recovered (about 2·0× longer in combined body and tail length). This raises a possibility that the diversification of paramphistomoids is more recent than the longer evolutionary history of Bulinus in Africa might suggest. More data are needed to resolve the phylogenetic position of this and other paramphistomoid clades, including those found in non-ruminant species.

In Kenya, Bulinus globosus, B. nasutus, B. africanus, B. tropicus, B. forskalii and Biomphalaria pfeifferi, are known to transmit paramphistomoids as well as ruminant and/or human schistosomes (Southgate et al. 1989; Brown, 1994). The overlap in use of snail hosts creates opportunities for distinctive interactions between the two common digenean groups. For example, in Kenya, Southgate et al. (1989) found that Bulinus tropicus was capable only of supporting the development of Schistosoma bovis to production of cercariae if it was first exposed to C. microbothrium. Similarly, in South America, Biomphalaria oligoza and Biomphalaria orbignyi are naturally resistant to S. mansoni, but become susceptible to S. mansoni if first exposed to Zygocotyle lunata (Spatz et al. 2012). Paramphistomoids can also have the opposite influence on the success of other digeneans during co-infections. For example, as compared to snails exposed only to F. hepatica, significantly fewer Pseudosuccinea columella produced F. hepatica cercariae if first exposed to Calicophoron daubneyi and then later exposed to F. hepatica (Dreyfuss et al. 2016).

This study has shown that even in a fairly circumscribed area within one East African country that a considerable diversity of paramphistomoid flukes is present and that several of these fluke species are abundantly represented. Paramphistomoids are of veterinary interest because of their ubiquitous presence in herds of cattle, sheep and goats that are routinely watered in natural habitats where the presence of susceptible species of snails ensures their transmission. Whether the species we have encountered have long parasitized domestic livestock or represent recent acquisitions from the region's many wild ruminants is an interesting question for future study. Studies currently underway in Kenya indicate that paramphistomoid infections are very common in some snail populations, so much so that they may represent significant impediments to the ongoing transmission of schistosomes using the very same snail hosts in the same aquatic habitats (Laidemitt M.R., personal communication, 2016). Furthermore, the spectra of freshwater snails used by these two common digenean groups are broadly overlapping, further increasing the likelihood that interesting interactions and accommodations have been made over evolutionary time. It will be interesting to more fully ascertain how these two major groups of digeneans influence one another's abundance. It is clear though that the domestication of livestock ensures that both paramphistomoid and schistosome (both human and ruminant schistosome species) life cycles are perpetuated side-by-side in the same habitats year after year. Livestock domestication may well prove to have had multiple downstream effects – mediated by the digeneans of livestock – on the present-day transmission of the all-too-common human blood flukes of sub-Saharan Africa.