COX1: Cytochrome oxidase 1
PCR: Polymerase chain reaction
ITS1: Internal transcribed spacer I
ITS2: Internal transcribed spacer II
mtDNA: Mitochondrial DNA
rDNA: ribosomal DNA
MBST: Molecular Biological System Transfer
Pulex irritans belongs to Siphonaptera as a relatively-small order of wingless holometabolous insects, known as ectoparasites in their adult stage [ 1 ].
Certain flea species can transmit pathogens such as Yersinia pestis, Rickettsia typhi, and Bartonella henselae through their bites, feces and saliva. Biting can also cause severe itching and skin infections [ 2 , 3 , 4 ]. P. irritans can cause serious health complications and transmit diseases in many parts of the world [ 5 ]. Global warming, international trade, travel, and population growth influence the epidemiology of these infections [ 6 , 7 , 8 ].
Cladistic analysis using morphological features including head, chest, and abdomen have been conducted to classify various fleas, nevertheless, accurately identifying certain flea species is difficult owing to variations in their morphological characteristics. Molecular markers are therefore used to accurately identify different flea species and investigate their phylogenetic relationships [ 9 , 10 , 11 , 12 ].
Molecular markers can be used to effectively characterize and specify an insect owing to their high stability and conserved nature; for instance, ribosomal DNA (rDNA) with its subunits such as 28s, 18s, 5.8s, ITS1 and ITS2 is commonly used to identify insect species [ 9 , 10 , 11 ]. As the largest protein-encoding gene in the mitochondria of metazoan organisms, the cytochrome oxidase subunit I (COXI) gene has been extensively used to perform phylogenetic research and identify species and their differences [ 13 , 14 ].
Accurately identifying flea species is crucial given that they constitute a vector of dangerous pathogens and thus a serious threat to public health. This study employed ITS1, ITS2, and COX1 as molecular markers to investigate the molecular characteristics of different populations of P. irritans and their phylogenetic relationships in four areas in Iran.
After morphological examinations using reliable identification keys for Iranian fleas, 1019 (52.61%) out of 1937 samples were identified as P. irritans (Table 1). Furthermore, the amplicons obtained from PCR for ITS1, ITS2 (rDNA), and COX1 (mtDNA) were approximately 1000 bp, 500 bp and 700 bp, respectively. The Clustal Omega based comparison of the ITS1 sequence in P. irritans showed a similarity of 99.59% and a molecular diversity of 0.41% between the sequences in ten different regions (Figure 1). A 100% similarity was observed among the sequences of Sanandaj (MN684797.1), Khorramabad (MN684787.1), Urmia (MN684784.1), Bahar (MN684777.1), and Kermanshah(MN 684776.1). These sequences were also similar (99.79%) to those of Hamadan (MN684779.1), Kuhdasht (MN684790.1), Kamyaran (MN684786.1), Gilangharb(MN684775.1), and Mahabad (MN684778.1) with four nucleotide differences. Four nucleotide polymorphisms were observed at positions 152, 165, 193 and 279 within the ITS1 spacer in the ten populations. Sequencing the ITS1 fragment showed three 99-bp replicas in the first half of the sequence beginning at positions 146, 246, and 352 and ending at position 431. Comparing the ITS2 gene sequences in the ten regions showed a 100% similarity between the sequences (Figure 2). Comparing the COX1 sequence of the nuclear genome (mtDNA) in the four different geographical zones (accession ID: MN173748-761) found all the species to belong to P. irritans and showed a 99.86% similarity between the sequences (Figure 3).
|Isolation source||Location||Longitude, latitude||Total||host||Accession number|
|Kermanshah||Ker||34.1397° N, 45.9206° E||86(8/4%)||41(47/6)||19(22)||4(4/6)||16(18/6)||6(6/9)||MN173752.1||MN684755.1||MN684792.1|
|Ghi||34.3277° N, 47.0778° E||84(8/4%)||30(35/7)||17(20/2)||15(17/8)||7(8/3)||15(17/8)||MN173749.1||MN684776.1||MN684782.1|
|Kurdestan||San||35.3219° N, 46.9862° E||92(9%)||22(23/9)||25(27/1)||20(21/7)||21(22/8)||4(4/3)||MN173760.1||MN68797.1||MN684801.1|
|Kam||34.7956° N, 46.9368° E||62(6%)||26(41/9)||3(4/8)||9(14/5)||24(38/7)||-||MN173754.1||MN684796.1||MN684781.1|
|West Azerbaijan||Urm||37.5498° N, 45.0786° E||139(13/64%)||58(41/7)||45(32/3)||6(4/3)||30(21/5)||-||MN173761.1||MN684784.1||MN684780.1|
|Mah||36.7684° N, 45.7337° E||63(6/1%)||4(6/3)||30(47/6)||13(20/6)||14(22/2)||2(3/1)||MN173758.1||MN684778.1||MN684803.1|
|Hamadan||Ham||34.9083° N, 48.4393° E||79(7/75%)||15(18/9)||35(44/3)||3(3/7)||26(32/9)||-||MN173753.1||MN684779.1||MN684800.1|
|Bah||34.9083° N, 48.4393° E||101(9/91%)||19(18/8)||50(49/5)||2(1/9)||29(28/7)||1(0/9)||MN173751.1||MN684777.1||MN684798.1|
|Lorestan||Kho||33.4647° N, 48.3390° E||98(9/61%)||8(8/1)||33(33/6)||8(8/1)||44(44/8)||5(5/1)||MN173756.1||MN684787.1||MN684794.1|
|Kuh||33.5275° N, 47.6111°E||93(9/12%)||2(2/1)||45(48/3)||3(3/2)||40(43)||3(3/2)||MN173757.1||MN684790.1||MN684793.1|
|Khozestan||Ahv||31.3183° N, 48.6706° E||21(2%)||-||5(23/8)||10(47/6)||3(14/2)||3(14/2)||MN173750.1||-||-|
|Mazandaran||Sar||36.5659° N, 53.0586° E||31(3%)||-||15(48/3)||9(29)||2(6/4)||5(16/1)||MN173759.1||-||-|
|Kerman||Kerm||30.2839° N, 57.0834° E||15(1/47%)||3(20)||-||2(13/3)||5(33/3)||10(66/6)||MN173748.1||-||-|
|North Khorasan||Esf||37.4710° N, 57.1013° E||55(5/39%)||-||-||40(72/7)||-||15(27/2)||MN173755.1||-||-|
|Sampling sites/regions: West (Ker=Kermanshah, Ghi=Ghilangharb, Kam=Kamiyaran, San=Sanandaj, Ham=Hamedan, Bah=Bahar, Kuh= Kuhdasht, Kho= Khorramabad), South (Ahv= Ahvaz, Kem=Kerman (, North (Sar=Sari)North East (Esf= Esfarayen )and Northwest (Urm=Urmia, Mah= Mahabad,) of Iran.|
Multiple alignments of amino acids in COX1 showed a 99.86% similarity and only one difference at position 54 of the methionine isolate and other isoleucine isolates in Hamadan (Figure 4). Two subclades of the phylogenetic tree based on the similarities between COX1 sequences in the present study and registered sequences in GenBank comprised subclade A1, including P. irritans of the present and previous studies conducted in Iran and that isolated in China, Turkey, Jerusalem, Spain and Croatia. Subclade A2 included Ctenocephalides felis in Iran and C. canis in Turkey. Subclades A3 and A4 included C. orientis and C. canis, respectively. Iranian subclade B also included Nosopsyllus fasciatus and Xenopsylla cheopis (Figure 5).
Three subclades of the phylogenetic tree based on the similarities between ITS1 sequences in the present study and registered sequences in GenBank comprised subclade A, including the sequences of the present and previous studies conducted in Iran. Subclades B and C included C. canis and N. fasciatus, respectively (Figure 6).
Four subclades of the phylogenetic tree based on the similarities between ITS2 sequences in the present study and registered sequences in GenBank comprised subclade C, including the sequences of the present study. Subclades A and B included C. canis and X. cheopis, respectively . Subclade D also included N. fasciatus (Figure 7).
The morphological characteristics of P. irritans were consistent with the molecular findings. Analyzing the partial COX1 gene can help identify P. irritans species and assess their intraspecific similarity. ITS1 and ITS2 sequences could also be used as useful markers for species-level differentiation and intraspecific similarity determination. The COX1 gene is, however, more efficient than the ITS1 and ITS2 genes in detecting genus, species and intra-species similarities. The diversity at the genus level and even between members of the same genus caused by high levels of replacement in the Internal transcribed spacer counteracted their beneficial effects. Given this diversity, the sequences of the heterogenic versions of ITS fragments are determined in cloning plasmids and then specifically, which requires laborious and costly studies[ 15 ].
Molecular methods have enabled the identification of different ectoparasite species with high morphological similarities [ 16 ]. The present study employed phylogenetic and molecular approaches to compare P. irritans in four different geographical zones of Iran. The morphological characteristics of P. irritans were found to be consistent with those obtained in previous studies conducted in Iran to determine the flea fauna in different hosts and cities. [ 17 , 4 , 18 ].
The host specificity can affect the intraspecific genetic diversity given the high levels of intraspecific genetic variations in general parasite species, which cause their infestation of broad ranges of hosts [ 19 ]. Research suggests no specific hosts for fleas, as they prefer diverse types of the host, which demonstrates intraspecific genetic diversity. Hornok et al. (2018) reported diversity in synanthropic flea species such as C. canis and P. irritans by investigating mitochondrial sequences [ 20 ]. The present morphological data showed no differences between the P. irritans specimens collected from different geographical zones. The results of the morphological investigation were also consistent with those of the molecular exploration. The present study reported a 100% similarity between the nucleotide sequence of COX 1 gene in the P. irritans and the sequence of this gene in the P. irritans in similar studies conducted in Iran. [MF380389.1, MF380390.1, MF380391.1 (NCBI)]. Similarities to the sequences in Spain (LT797470.1) and China (MF000666.1) were also reported as 99.85% and 99.55%, respectively. Furthermore, the present findings were consistent with those obtained by Seyyedzadeh et al. (2018), who reported no differences between the isolates of C. canis in different areas of West Azerbaijan, Iran [ 21 ] and with those obtained by Zurita et al. (2019), who observed no significant differences in morphological data between P. irritans in Spain and Argentina. In addition, Zurita et al. (2019) reported a significant intraspecific similarity between the two populations based on mitochondrial genes [ 22 ]. Hornok et al. (2018) observed no morphological differences between human and wild carnivorous P. irritans specimens in Hungary and Croatia [ 20 ], which is consistent with the present results. In contrast, Krasnov et al. (2015) found morphological differences between the flea species isolated from different hosts in different geographic areas to suggest high levels of genetic diversity [ 23 ].
Evidence suggests ITS1 and ITS2 constitute appropriate molecular markers for analyzing phylogenetic relationships in fleas at their species level [ 24 , 11 ]. The present research observed a 99.59% similarity between the nucleotide sequences of the ITS1 spacer, in five provinces in Iran. Analyzing the nucleotide sequence of the ITS1 spacer, also showed three 99-bp replicas in its first half, which is consistent with the results of Ghawami et al. (2018), who compared the nucleotide sequence of the ITS1 spacer, in P. irritans in two geographic areas and reported only one nucleotide difference with a 99.85% similarity [ 25 ]. In line with the present results, those obtained by Vobis et al. (2004) suggested a relatively-constant nucleotide sequence of the ITS1 spacer, in different C. felis populations [ 11 ].
A 100% similarity was observed between the ITS2 nucleotide sequence of different samples. Zurita et al. (2019) reported an intraspecific variation of 90.1%-100% in four C. canis clones based on the ITS2 spacer,, which is consistent with the present results. Vobis et al. (2004) also used the ITS2 sequence to identify different flea species and determine their intra-species differences [ 11 ]. P. irritans and C. canis belonging to the Pulicidae family constitute the main flea genera infesting human and cattle in the west and northwest of Iran. No significant differences were observed among the samples collected from different hosts. Insignificant differences were found between the study populations in terms of mitochondrial markers and the nuclear ribosome. The results of plotting the phylogenetic tree showed all the samples collected from different regions in the same branch despite their negligible sequence differences. Analyzing mitochondrial genome also showed that this sequence is more useful than the nuclear ribosomal genome in illustrating intra-species similarity and differentiation at genus and species levels.
The Iranian plateau is located in the northern hemisphere at 32.4279 °N, 53.6880 °E. The four study geographical zones in Iran included region 1: Caspian Sea (temperature: 8.00-26.00 ˚C, annual precipitation: 400-1500 mm), region 2: Mountainous area (temperature: -5.00-29.00 ˚C, annual precipitation: 200-500 mm), region 3: Persian Gulf (temperature: 12.60-35.00 ˚C, annual precipitation: 200-300 mm) and region 4: The Central Desert (temperature: -4.00-44.00 ˚C, annual rainfall: below 100 mm). According to meteorological data, Iran is climatically divided into three zones, i.e. warm-dry, cold, temperate-humid, and warm-humid [ 26 ]. Four different geographical zones in this country from which adult fleas were collected included Kermanshah, Kordestan, Hamadan, Lorestan, West Azerbaijan, Khuzestan, Khorasan, Mazandaran, and Kerman (Figure 8).
A total of 1937 flea samples were directly isolated from the host body and houses using light traps and human-baited traps. They were transferred in 70% ethanol to the parasitology laboratory of the Faculty of Veterinary Medicine. They were then cleaned with distilled water several times, immersed in 5% potassium hydroxide for 24 hours, re-washed with distilled water, and dehydrated with a graded series of alcohol (50%-96%(. The fleas were ultimately identified by observing the slides with an optical microscope according to diagnostic keys [ 27 , 28 ].
All the P. irritans samples collected from different provinces and hosts were morphologically characterized by the following features: the round-headed anterior margin with a pair of eyes and no protrusion, absence of genal and pronotal ctenidia, asymmetrical club of antenna, strong setae under the eyes, row of small spines mass on the inner surface of the posterior side of the coxa in females 7-10 and in males 8-12, and no pleural rod on the mesopleural segment [ 27 ]. Comparing the samples collected from different provinces and hosts based on morphological characteristics identified 52.61% to belong to P. irritans species despite their negligible differences.
The samples were individually crushed with a sterile scalpel, poured into a sterile microtube and their total DNA was extracted using DNA extraction kits (MBST, Tehran, Iran) according to the manufacturer’s instructions. The extracted DNA samples were stored in sterile microtubes at -20 °C until use. DNA was extracted from the samples isolated from different hosts in different areas. The specific primers used by Vobis et al (2004). were employed to identify P. irritans and 1000-bp ITS1 and 500-bp ITS2 fragments in the rDNA of P. irritans. Primers used by Folmer et al. (1994) were used to amplify the mitochondrial COX1 gene fleas P. irritans (Table 2). The result of this proliferation was observed as a 700 bp band [ 29 ]. To determine molecular differences between the species isolated from different hosts and areas, PCR was performed in a total volume of 50 µl, containing a 5-µl DNA template, a5-µl 10x PCR buffer, 1-µl Taq Polymerase, 1 µl of each primer (20 µM), 1 µl of dNTP (100 µM) and 4 µl of MgCl2 (50 mM). All these materials were procured from SinaClon, Iran. The thermal cycle of the PCR on fragments ITS1 and ITS2 was as follows: incubation at 94 °C for 5 minutes to denature double-stranded DNA followed by 30 cycles of one-minute denaturation cycles at 94 °C, annealing at 54 °C for 1 minute, extension at 72 °C for 1 minute and an additional extension at 72 °C for 5 minutes. The COX1 gene also underwent replication at 95 °C for 5 minutes followed by 35 cycles of denaturation cycles at 95 °C for 1 minute, annealing at 55 °C for 45 seconds, extension at 72 °C for 45 seconds, and final extension at 72 °C for 10 minutes. The PCR products were electrophoresed on 1.5% agarose gel, stained with a safe stain, visualized under UV light and then purified for sequencing. After purifying the PCR products, 34 samples, including 10 positive samples for ITS1, 10 positive for ITS2, and 14 positive for COX1 sequence, collected from different hosts were transferred to Takapouzist Company for sequencing.
|Gene||Primer nucleotide sequence||Annealing temperature (°C)||Fragment length (bp)||Reference|
|ITS1||Sen: GTA CAC ACC GCC CGT GCG TAC T||54 °C||1000 bp||Lawrence et al. 2014|
|Rev: GCT GCG TTC TTC ATC GAC CC|
|ITS2||Sen: GGG TCG ATG AAG AAC GCA GC||54 °||500 bp||Lawrence et al. 2014|
|Rev: GCG CAC ATG CTA GAC TCC GTGGTT CAA G|
|COX1||Sen: GGT CAA CAA ATC ATA AAGATA TTG G||55 °C||700 bp||Folmer et al.1994|
|Rev: GAA GGG TCA AAG AAT GAT GT|
A total of 34 PCR products were transferred to Takapouzist for sequencing. An accession number was obtained by online recording the sequencing results on NCBI by the host and geographical area. The data were analyzed using blasting sequences on NCBI and plotting the phylogenetic tree in Mega 6 using maximum likelihood and bootstrapping (1000 replicates). The clustering method proposed by Zurita et al.(2019) was also conducted. The required sequences were ultimately aligned and compared using EMBOSS Needle and Clustal Omega to evaluate the similarity percentage of the nucleotide sequences of mitochondrial (mtDNA) and ribosomal genomes (rDNA) in different provinces and hosts.
Conceived and designed the experiments and revised the manuscript draft: S.S., M.T., performed the experiments, analysed the data and drafted the manuscript: S.S. All authors approved the final version of the manuscript.
The authors would like to express their gratitude to the authorities of the Faculty of Veterinary Medicine and Research Council of Urmia University for their approval and financial support of this research.
The authors declare that there is no conflict of interest.
- Lewis RE. Resume of the Siphonaptera (Insecta) of the world. J Med Entom. 1998; 1; 35(4): 377-389.
- Azizi MH, Azizi F. A history of the human plague in Iran. Arch Iran Med. 2010; 13(6): 563-569.
- Esamaeili S, Azadmanesh K, Naddaf SR, Rajerison M, Carniel E, Mostafavi E. Serologic survey of plague in animals, Western Iran. Emerg Infec Dis. 2013; 19(9):1549.
- Tavassoli M, Ahmadi A, Imani A, Ahmadiara E, Javadi S, Hadian M. Survey of flea infestation in dogs in different geographical regions of Iran. Korean J parasitol. 2010; 48(2):145.
- Bitam I, Dittmar K, Parola P, Whiting MF, Raoult D. Fleas and flea-borne diseases. In J Infec Dis. 2010; 14(8): 667-676.
- Leulmi H, Socolovschi C, Laudisoit A, Houemenou G, Davoust B, Bitam I, Raoult D, Parola P. Detection of Rickettsia felis, Rickettsia typhi, Bartonella species and Yersinia pestis in fleas (Siphonaptera) from Africa. PLoS Negl Trop Dis. 2014; 8(10): 3152.
- Ratovonjato J, Rajerison M, Rahelinirina S, Boyer S. Yersinia pestis in Pulex irritans fleas during plague outbreak, Madagascar. Emerg Infect Dis. 2014; 20(8):e1414.
- Stenseth NC, Atshabar BB, Begon M, Belmain SR, Bertherat E, Carniel E, Gage KL, Leirs H, Rahalison L. Plague: past, present, and future. PLoS Med. 2008; 5(1):e3152.
- Azarm A, Dalimi A, Mohebali M, Mohammadiha A, Zarei Z. Morphological and molecular characterization of Ctenocephalides spp isolated from dogs in north of Iran. J Entomol. 2016; 4(4):713-717.
- Gamerschlag S, Mehlhorn H, Heukelbach J, Feldmeier H, D’Haese J. Repetitive sequences in the ITS1 region of the ribosomal DNA of Tunga penetrans and other flea species (Insecta, Siphonaptera). Parasito Res. 2008; 102(2): 193-199.
- Vobis M, D’Haese J, Mehlhorn H, Mencke N, Blagburn BL, Bond R, Denholm I, Dryden MW, Payne P, Rust MK, Schroeder I. Molecular phylogeny of isolates of Ctenocephalides felis and related species based on analysis of ITS1, ITS2 and mitochondrial 16S rDNA sequences and random binding primers. Parasitol Res. 2004; 94(3): 219-226.
- Zurita A, Callejon R, Garcia‐sanchez ÁM, Urdapilleta M, Lareschi M, Cutillas C. Origin, evolution, phylogeny and taxonomy of Pulex irritans. Med Vet Entomol. 2019; 33(2): 296-311.
- Bowen BW, Clark AM, Abreu-Grobois FA, Chaves A, Reichart HA, Ferl RJ. Global phylogeography of the ridley sea turtles (Lepidochelys spp.) as inferred from mitochondrial DNA sequences. Genetica. 1997; 101(3): 179-189.
- Simon C, Frati F, Beckenbach A, Crespi B, Liu H, Flook P. Evolution, weighting, and phylogenetic utility of mitochondrial gene sequences and a compilation of conserved polymerase chain reaction primers. Ann Entomol Soc Am. 1994; 87(6): 651-701.
- Mclain DK, Wesson DM, Oliver Jr JH, Collins FH. Variation in ribosomal DNA internal transcribed spacers 1 among eastern populations of Ixodes scapularis (Acari: Ixodidae). J Med Entomol. 1995; 32(3):353-560.
- Lv J, Wu S, Zhang Y, Chen Y, Feng C, Yuan X, Jia G, Deng J, Wang C, Wang Q, Mei L. Assessment of four DNA fragments (COI, 16S rDNA, ITS2, 12S rDNA) for species identification of the Ixodida (Acari: Ixodida). Parasites vectors. 2014 Mar 3; 7:93.
- Bahrami S, Kharrati L, Makki M. The First Case Report of Human Infestation with Ctenocephalides canis in Behbahan (Khuzestan Province). J of Ardabil Univ of Med Sci. 2016; 16(3): 317-322.
- Yakhchali M, Hosseine A. Prevalence and ectoparasites fauna of sheep and goats flocks in Urmia suburb, Iran. Vet Arhiv. 2006; 76(5): 431-442.
- Van der Mescht L, Matthee S, Matthee CA. Comparative phylogeography between two generalist flea species reveal a complex interaction between parasite life history and host vicariance: parasite-host association matters. BMC Evol Biol. 2015; 15(1):1-5.
- Hornok S, Beck R, Farkas R, Grima A, Otranto D, Kontschan J, Takacs N, Horváth G, Szoke K, Szekeres S, Majoros G. High mitochondrial sequence divergence in synanthropic flea species (Insecta: Siphonaptera) from Europe and the Mediterranean. Parasites vectors. 2018; 11(1):1-1.
- Seyyed‐Zadeh SJ, Bozorg‐Omid F, Telmadarraiy Z, Terenius O, Chavshin AR. Evidence for the presence of Ctenocephalides orientis in livestock dwellings in northwest Iran. Med Vet Entomol. 2018; 32(3): 383-387.
- Zurita A, Callejon R, De Rojas M, Halajian A, Cutillas C. Ctenocephalides felis and Ctenocephalides canis: introgressive hybridization?. Syst Entomol. 2016; 41(3): 567-579.
- Krasnov BR, Shenbrot GI, Khokhlova IS. Historical biogeography of fleas: the former Bering Land Bridge and phylogenetic dissimilarity between the Nearctic and Palearctic assemblages. Parasitol Res. 2015; 114(5): 1677-1686.
- Marrugal A, Callejón R, de Rojas M, Halajian A, Cutillas C. Morphological, biometrical, and molecular characterization of Ctenocephalides felis and Ctenocephalides canis isolated from dogs from different geographical regions. Parasitol Res. 2013; 112(6): 2289-2298.
- Ghavami MB, Mirzadeh H, Mohammadi J, Fazaeli A. Molecular survey of ITS1 spacer and Rickettsia infection in human flea, Pulex irritans. Parasitol Res. 2018; 117(5):1433-1442.
- Skerman KD, Shahlapoor AA, Eslami AH, Eliazian M. Observations on the incidence, epidemiology, control and economic importance of gastro-intestinal parasites of sheep and goats in Iran. Arch Razi Inst. 1970; 22(1): 187-196.
- Asmar M, Piazak N, Karimi Y. An Illustrated Key for Flea of Iran(Persian). Pasteur Institute of Iran, Research note. 1979;2-15.
- Shayan A, Rafinejad J. Arthropod parasites of rodents in Khorram Abbad district, Lorestan Provincen of Iran. Iranian J Public Health. 2006; 70-6.
- Lawrence AL, Hii SF, Jirsova D, Panakova L. Integrated morphological and molecular identification of cat fleas (Ctenocephalides felis) and dog fleas (Ctenocephalides canis) vectoring Rickettsia felisin central Europe. Vet Parasitol. 2015; 10(3-4): 215-23.
- Folmer O, Black M, Hoeh W, Lutz R, Vrijenhoeket R. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Mol Ma Biol Biotechnol. 1994; 3(5): 294-299.