Molecular characterization and phylogeny of Linguatula serrata (Pentastomida: Linguatulidae) based on the nuclear 18S rDNA and mitochondrial cytochrome c oxidase I gene.

Linguatula serrata, a cosmopolitan parasite, is commonly known as tongue worm belonging to the subclass Pentastomida.We collected the nymphal stage of the worm from mesenteric lymph nodes of cattle and identified these as L. serrata based on morphology and morphometry. The 18S rDNA sequences showed no intraspecific variation, although cox1 sequences showed 99.7-99.9% homology. In the phylogenies inferred from both gene loci, members of the genus Linguatula (order Porocephalida) were closer to those of the order Cephalobaenida than to those of Porocephalida, reflecting a mismatch with the corresponding morphology-based taxonomy. Accordingly, analyses of additional gene loci using a larger number of taxa across the Pentastomida should be undertaken to determine an accurate phylogenetic position within the Arthropoda.

Linguatula serrata Fröhlich, 1789, belonging to the subclass Pentastomida (tongue worms), typically inhabits in the upper respiratory tract (nasal airways and sinuses) of carnivorous mammal [6, 13] and has a worldwide distribution. The worm has an indirect life cycle where a wide range of mammals, including cattle, goat, sheep, camels, and rodents and rabbits act as intermediate hosts. The infective stage, nymph, usually occurs in the mesenteric lymph nodes, liver and lungs of intermediate hosts, and the final hosts become infected through ingestion of raw or improperly cooked meat containing the nymph, or by the drinking water contaminated with the feces or nasal excreta of carnivores containing eggs. Humans harbor both the adult [3, 13] and nymphal stages of the worm, which can cause nasopharyngeal and visceral linguatuliasis, respectively [18]. Ocular linguatuliasis, a very rare form, in human has also been reported from the United States [2], Congo [14], Ecuador [8] and Austria [7] with the consequences ranging from photophobia and exophthalmia to permanent blindness. Therefore, L. serrata bears public health significance worldwide with severe consequences, especially in the developing world. Although there have been many reports on the incidence and morphology [6, 8, 11, 12, 15] and diagnosis [18] of L. serrata, there have been only two reports on the molecular properties of the worm [5, 7], and the phylogenetic relationship of pentastomids has yet to be resolved. The nuclear 18 subunit ribosomal DNA (18S rDNA) has proved to be a potential marker in the barcoding approach to characterizing eukaryotes [4, 9] and copepods [22]. Its highly conserved nature and slow rate of evolutionary mutation make 18S rDNA suitable for inter-species distinction at higher taxonomic levels [17]. In contrast, the mitochondrial cytochrome c oxidase subunit 1 (cox1) gene has recently been used for phylogenetic and population genetics analyses in different insect species based on its high mutational rate compared to 18S rDNA. In this study, we analyzed the nucleotide sequences of 18S rDNA and cox1 of L. serrata, and attempted to elucidate the molecular phylogenetic relationships of this species within the Pentastomida and with other taxonomically related organisms.

In January 2015, a total of nine pentastomid nymphs were recovered from the mesenteric lymph nodes of seven cattle at slaughterhouses in Bogra, Khulna and Mymensingh in Bangladesh. The nymphs were fixed in 70% ethanol. Two of the specimens were immediately stained with hematoxylin-carmine solution to enable photography of the whole organisms, and these were not included in the subsequent molecular studies. The remaining seven nymphs were measured for the body length, maximum body width (MW) and distance from the anterior extremity to the MW under an optical microscope. Tissue from the left side of the intestine of each specimen was then removed for DNA extraction. The remaining body parts were stained with hematoxylin−carmine solution for detailed morphometric and morphological studies.

Total DNA was extracted from the nymphs using a High Pure PCR template preparation kit (Roche, Mannheim, Germany) following the instructions of the manufacturer, and stored at −20°C for future use. DNA fragments of 18S rDNA and cox1 were amplified by the polymerase chain reaction (PCR) in a total volume of 10 µl, containing 0.25 U TKs Gflex™ DNA polymerase (TaKaRa Bio Inc., Otsu, Japan), 0.2 µl template DNA, 1 µM each primer and 2× Gflex PCR buffer. The primer sets used were SSU1 [21] and P2rev [20] for the 18S rDNA fragment (1,000 bp) and Cox1 LF and Cox1LR [5] for cox1 (1,015 bp). The reaction cycle for 18S rDNA consisted of an initial denaturation at 94°C for 1 min; followed by 40 cycles at 98°C for 10 sec, 60°C for 15 sec and 68°C for 30 sec, whereas that for cox1 consisted of an initial denaturation at 95°C for 15 min; followed by 45 cycles at 94°C for 45 sec, 47°C for 45 sec and 72°C for 90 sec; with a final extension at 72°C for 10 min. The PCR reactions were performed using a GenAmp PCR system 2700 (Applied Biosystems Japan Ltd., Tokyo, Japan). PCR products were separated by 1.2% agarose gel electrophoresis, stained with ethidium bromide and visualized with an ultraviolet transilluminator. The PCR products were subsequently purified by gel extraction using a NucleoSpin Gel and PCR clean-up kit (MACHEREY−NAGEL GmbH & Co. KG, Düren, Germany).

The PCR amplicons for both 18S rDNA and cox1 were sequenced directly from both directions using a BigDye Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems, Foster City, CA, U.S.A.). For sequencing, a pair of additional primers, P1fw and P1rev [20], were used in case of 18S rDNA, whereas UK-F (5′-CACTGCTTTCCTCCTAATCCTCTC-3′) and UK-R (5′-GTGTTGGAAGAGGATTGGGTC-3′) were used for cox1. The thermal program used for the sequencing reaction comprised an initial denaturation at 96°C, followed by 25 cycles at 96°C for 10 sec, 56°C for 5 sec and 60°C for 4 min. The resultant amplicons of both 18S rDNA and cox1 were then precipitated with ethylenediamine tetra acetic acid (EDTA)/ethanol for final purification prior to sequencing. Sequence reactions were run in a 3500 Genetic Analyzer (Applied Biosystems). The sequence data were initially assembled using ATGC ver. 6.0.3 (Genetyx Co., Tokyo, Japan) and were analyzed using GENETYX version 10.0.2 (Fluxus Technology Ltd., Suffolk, U.K.).

Haplotypes for both 18S rDNA (1000 bp) and cox1 (1015 bp) were distinguished by aligning the respective sequences in Clustal X version 2.0 [19]. The phylogenetic analyses were conducted based on the sequences of 18S rDNA (913 bp) and cox1 (605 bp) by Maximum Likelihood (ML) using MEGA6 [16]. On the basis of the Bayesian Information Criterion for ML analyses, the best nucleotide substitution model was determined. The Kimura 2 parameter model with discrete Gamma distribution (K2+G) with 5 rate categories was selected for constructing the ML tree inferred from the 18S rDNA sequences, whereas the Hasegawa–Kishino–Yano model with discrete Gamma distribution (HKY+G) was selected for the cox1 sequences. All the positions containing gaps or missing data were eliminated. For construction of the ML tree inferred from 18S rDNA sequences, 21 reference sequences were included: 15 pentastomid sequences, L. serrata (JX088397, KT581431, KT581433, KP100453 and FJ528908), L. arctica (KF029439 − 42), Armillifer agkistrodontis (FJ607339), A. moniliformis (HM048870), Raillietiella sp. (AY744887 and EU370434), Reighardia sternae (AY304521) and Hispania vulturis (AY304520); three crustacean taxa, Hepatus epheliticus (AF436004), Squilla empusa (L81946) and Caligus elongatus (AY627020); two horse-shoe crabs (Xiphosura), Limulus polyphemus (HQ588740) and Tachypleus gigas (HQ876477); and one insect taxon, Thermobia domestica (AF370790). A total of 17 references sequences were included for construction of the ML tree inferred from cox1: nine other pentastomid sequences, L. serrata (KF029447), L. arctica (KF029443 − 46), A. agkistrodontis (FJ607340), A. armillatus (AY456186, NC_005934), Raillietiella sp. (JF975594); three nematode taxa, Trichinella spiralis (AF293969), Trichuris ovis (JQ996232) and T. trichuria (GU385218); two crustacean taxa, Speleonectes tulumensis (AY456190) and Squilla empusa (DQ191684); two insect taxa, Ceratitis curvata (AY788423) and T. domestica (AF639935); and one trematode taxon, Fasciola gigantica (NC_024025). In both the trees, bootstrap analyses were conducted using 1,000 replicates.

The nymphs are white, dorso-ventrally flattened, vermiform with a broader and rounded anterior extremity, and tapering at the posterior end (Fig. 1A). The cuticle of the nymphs has on an average 88 rings (80 − 94), and each segment contains a row of spines (Fig. 1C and 1D) at its posterior border along the body. The length and number of spines (Fig. 2E) vary in the anterior, middle and posterior regions. The interval between the rings also varies in different areas of the body, with the largest interval being between the first and second rings (Fig. 1B and 1E; Table 1). The buccal framework, which is rectangular with rounded corners and chitinous borders, is located ventro-medially at the level of the first ring in the anterior region of the body. Two pairs of hooks (Fig. 1B) with sharp tips (blade) and an almost rectangular base (handle) are located in the anterior region and have an area of approximately 11,113.6 µm2. The first pair of hooks are located just lateral to the buccal framework at the level of the first ring, whereas the second pair are located postero-lateral to the buccal framework at the level of the second ring. In the region posterior to each hook, there is a rectangular and somewhat curved structure referred to as the fulcrum (Fig. 1A and 1B). Additionally, there are four accessory hooks anterior to the main hooks. A long slender intestine is located posterior to the buccal framework, and this extends to the posterior end along the median line, reaching to the terminal anus (Fig. 1A). Accompanying the intestine, a red primordial uterus extends from the anterior region to the posterior region in a coiling fashion. All of these morphological features and morphometric data (Table 1) are consistent with those of the previously described isolates of L. serrata [8, 12, 18], and therefore, the nymphs were morphologically identified as L. serrata.

Fig. 1.

Different body parts of a Linguatula serrata nymph. A: Representative entire nymph of Linguatula serrata (ventral view) from the lymph node of cattle. The broken lines represent the portion of the worm used for DNA extraction. Although this specimen was not included in the morphometric and molecular analyses, it was morphologically identical to those nymphs that were analyzed. B: Ventral view of the anterior region (head) of a Linguatula serrata nymph. “*” indicates accessory hooks. “**” indicates the main hooks. The white arrow indicates the buccal framework. The black arrow shows the method of measuring the interval between the first and second rings. The dotted lines designate the structural area of hooks and accessory hooks. The solid lines represent the points of measurement. C: Posterior part of the worm. D: showing the arrangement of ring spines along the posterior border of each segment/ring. E: showing the method of measuring the spine length (downward arrow) and ring interval (upward arrow).

Fig. 2.

Maximum likelihood trees for Linguatula serrata isolates recovered from cattle in comparison with the sequences of other pentastomids and arthropods deposited in GenBank. A: tree based on 18S rDNA sequences (913 bp). Bootstrap values higher than 80% are shown on the tree nodes. B: tree based on cox1 sequences (605 bp). Bootstrap values higher than 70% are shown on the tree nodes. Genotypes in bold are those generated in the present study.

Table 1.

Morphometric data of Linguatula serrata nymph

Region measured		This study			Lazo et al. [7]	Rezaei et al. [11]
		Minimum	Maximum	Average
Body (mm)	Length	3.0	5.2	4.2	4.7	4.9
	Maximum width (MW)	0.96	1.3	1.2	1.0	1.0
	Distance from anterior extremity to MW	1.0	1.4	1.3
Buccal framework	Length (mm)	0.16	0.18	0.17	0.19	0.15
	Width (mm)	0.069	0.11	0.089	0.092	0.079
	Area (µm2)	11,656.99	20,432.13	15,610.47	13,673.80
Hook (mm)	Blade	0.17	0.25	0.21	0.198
	Handle	0.090	0.15	0.13	0.14
Accessory hook	Length (mm)	0.069	0.14	0.10	0.096
	Width (mm)	0.020	0.047	0.030
	Area (µm2)	1,551.23	4,871.14	3,252.92	2,013.40
Falcrum	Length (mm)	0.28	0.39	0.33	0.35
	Width (mm)	0.058	0.11	0.088
	Area(µm2)	29,172.50	35,894.28	29,172.50
No. of segment/ring on the cuticle		80	94	88	95
Interval between the first and second ring (mm)		0.10	0.13	0.12	0.15
Ring interval (μm)	Anterior	63.34	88.8	73.43	60.4
	Middle	36.92	53.79	44.12	53.8
	Posterior	29.33	47.79	38.95	47.4
Length of ring spine (μm)	Anterior	33.25	41.9	37.24	32.5
	Middle	0.021	0.029	0.025	0.019
	Posterior	0.017	0.024	0.022	0.016
No. of Spine/segment	Anterior	85	106	96.14
	Middle	80	112	92.85
	Posterior	32	44	36

The nucleotide sequences (1,000 bp) of 18S rDNA from all seven L. serrata showed complete identity and were designated as genotype LS_Bd (LC150568). The sequence of LS_Bd is identical to the 18S rDNA sequence of L. serrata from Norway (JX088397). The sequences of L. serrata (LC250568 and JX088397) had 99.9% identity to that of L. arctica (KF029439-42) and formed a monophyletic clade in the ML tree based on the 18S rDNA sequences. In contrast, LS_Bd showed 98.3% identity (with 15 nucleotide substitutions) to L. serrata (FJ528908), which was isolated from a human eye in Austria. These results indicated that L. serrata (FJ528908) either was misidentified morphologically or that the sequence generated had ambiguity, as mentioned in a previous study [5]. In the ML tree of 18S rDNA, the members (H. vulturis, R. sternae and Raillitiella sp.) of the family Cephalobaenidae of Pentastomida formed a sister group to the clade formed by members of the family Linguatulidae (Fig. 2A). All of the aforementioned pentastomids formed a sister clade (with maximum support) to the cluster formed by A. agkistrondontis (FJ607339) and A. moniliformis (HM048870) of the family Porocephalidae.

The nucleotide sequences (1015 bp) of cox1 displayed 99.7−99.9% identity with three nucleotide substitutions and yielded four genotypes designated LS_COI_1 to LS_COI_4 (LC150781−LC150784). Among the four genotypes, LS_COI_3 had only one nucleotide substitution compared with the cox1 sequence of L. serrata (KF029447). In the ML tree, these genotypes formed a monophyletic clade along with L. serrata (KF029447), which was sister to the L. arctica clade with maximum support value (Fig. 2B). Similar to the 18S rDNA tree, Linguatula (Porocephalida) in the cox1 tree is closer to Raillitiella sp. (Cephalobaenida) than to Armillifer, a member of the order Porocephalida. The pentastomids are classified into two orders, Cephalobaenida and Porocephalida. The genera, Linguatula and Armillifer, are included into the order Porocephalida, and Raillietiella, Hispania and Reighardia into Cephalobaenida [10]. This morphology-based taxonomy is inconsistent with the present molecular phylogenies, which are supported by previous observations [5, 7]. If the cox1 sequences of nematodes, crustaceans, C. curvata and F. gigantica are removed from the data set, and T. domestica is used as an outgroup for constructing the ML tree, Linguatula becomes closer to Armillifer, matching the morphology-based taxonomy (Fig. 3). Additionally, the phylogeny based on the cox1 sequences indicates that pentastomids are more closely related to Nematoda than to Arthropoda (Fig. 2B). This observation is supported by previous findings [1], which suggested that Pentastomida is an intermediary group between the Arthropoda and some Nemathelminthes. Unfortunately, there are insufficient molecular data for most of the pentastomids to enable a resolution of this issue. Consequently, revision on the morphology-based taxonomy or inclusion of molecular data from an extensive range of taxa will be necessary to elucidate the phylogenetic relationships of the Pentastomida.

Fig. 3.

Maximum likelihood tree inferred from cox1 sequences (605 bp) of Linguatula serrata isolates recovered from cattle in comparison with the sequences of other pentastomids deposited in GenBank. Genotypes in bold are those designated in the present study. Bootstrap values higher than 90% are shown on the tree nodes.