Mitochondrial DNA discriminates distinct population of two deadly snakes (Reptilia: Elapidae) in Northeast India

Abstract

The DNA data of Indian snakes are scanty in the global database, especially from the northeastern states. The present study generated the mitochondrial Cytochrome b gene information of two morphologically identified deadly elapid species from Mizoram. Both, the King Cobra (Ophiophagus hannah) and Banded Krait (Bungarus fasciatus) showed monophyletic clades in the BA topology and cohesively clustered with the database sequences generated from distant geographical locations. The studied O. hannah depicted 2.7–7.6% K2P genetic distances with the specimens collected from China, Vietnam, and Thailand. Further, the northeast Indian B. fasciatus revealed 3.3–4% K2P genetic distance from Chinese, Vietnamese, Thailand, Indonesian, and Australian specimens. The TCS network showed distinct haplotypes for both the species collected from northeast India. The genetic information of these venomous snakes would be helpful for further rapid identification from the museum as well as from road-killed specimens, curbing the venom poaching and medical avenues.

Keywords:

• Elapidae
• venomous
• mitochondrial DNA
• phylogeny
• haplotypes

Introduction

The legless reptiles, snakes are distributed in every continent except Antarctica with more than 3000 species till date (Wallach et al. 2014; Uetz and Hošek 2019). Among them, over 300 species are known from the Indian subcontinent of which 61 are venomous (Whitaker and Martin 2015). The venomous snakebites are a significant risk to human health around the world and most impacted to the rural peoples (Warrell 2010; WHO 2019). Due to the evolution of modified teeth and venom glands, these species are skilled to inject venom into the prey animals (Gong et al. 2010). Owing to the high mortality, disability, and priority research, the World Health Organization (WHO) enlisted the venomous snakebite as a neglected tropical disease (Warrell and WHO 2009). The probable number of snakebites was 5.8 million worldwide with 1.25 lakh deaths, of which over one million snakebites and 50,000 deaths occurred every year in India (Mohapatra et al. 2011; Dandona et al. 2018). In the recent past, WHO launched the snakebite envenoming guidelines with the aspiration of diminishing the number of deaths and disabilities by 50% before 2030 (Bolon et al. 2019; WHO 2019). In a medical point of view, most of the snakes are nonvenomous; however, the deadly snakes are categorized into two groups: the highly noxious, common, widespread species, and the species with restricted distribution remote from human populations (WHO 2010; Williams et al. 2019).

In India, three Elapidae species (Bungarus caeruleus, Naja kaouthia, Naja naja) and three Viperidae species (Daboia russelii, Echis carinatus, Hypnale hypnale) are grouped into first category; however, seven Elapids (Bungarus fasciatus, Bungarus niger, Bungarus sindanus, Bungarus walli, Naja oxiana, Naja sagittifera, Ophiophagus hannah) and five Vipers (Cryptelytrops albolabris, Cryptelytrops purpureomaculatus, Trimeresurus malabaricus, Trimeresurus gramineus, Macrovipera lebetina) are grouped into second category (Alirol et al. 2010; Menon et al. 2017). Nevertheless, without proper identification of species and lack of awareness, the majority of victims are immediately scared, increasing heart beat and feeling faint before medical intervention. Nevertheless, 80% of the world’s population including Indians rely on traditional medicine for snakebite envenoming and primary health-care. The rural people of India largely believe in medicinal plant-based traditional practices, juice or paste of leaves, stem barks, and roots that are applied externally or consumed orally for snakebite first aid (Goswami et al. 2014; Lalramnghinglova 2016). Apart from this, the accurate identification of biting species is crucial for appropriate treatment, especially for neurotoxic envenoming. It is often encountered that the victim does not respond to the antivenom prescribed by the medical person due to the lack of species identification (Bhattacharya and Chakraborty 2007). Additionally, the manufacture of antivenoms and their effectiveness against the venoms of geographical widespread snakes is still underway (Ghosh et al. 2016). Hence, the knowledge of indigenous species, their genetic structure, the chemistry of venom, and signs, as well as symptoms of envenomation, is crucial.

The DNA-based species identification is largely used as an effective tool for systematics research, biodiversity assessment, and forensic sciences (Hebert et al. 2003; Dubey et al. 2011; Tyagi et al. 2019). These molecular techniques are tremendously used to identify the snakes from widespread geographical regions (Nagy et al. 2012; Chambers and Hebert 2016). Further, the mitochondrial DNA is also evident to identify the snake species from venom samples (Singh et al. 2012; Supikamolseni et al. 2015; Sharma et al. 2016; Smith et al. 2018). Besides, the phylogeny, population structure, and evolutionary history of snakes were effectively adjudicated by the assessment of both mitochondrial and nuclear genes (Slowinski and Lawson 2002; Pyron et al. 2013; Figueroa et al. 2016). Hence, in the present study, we sampled two deadly Elapidae snakes from Mizoram state in northeast India. We used both morphology and mitochondrial Cytochrome b sequence information to confirm the species identity. Further, we estimated the phylogenetic tree, genetic divergence, and haplotype to determine the status of these two snakes in northeast India. The present data would further be helpful to detect the road-kill species, arresting the venom smugglers racket, and medical avenues for specific antivenom production.

Materials and methods

The field survey and sampling of King Cobra, O. hannah (23.73 N 92.67E) (Figure 1A) and Banded Krait, B. fasciatus (24.01 N 92.37E) (Figure 1B) was conducted by the herpetology team of Mizoram University after obtaining the permission from Chief Wildlife Warden of Environment, Forests and Climate Change, Government of Mizoram. Both the specimens were morphologically identified and vouchered at the Department of Zoology, Mizoram University, and the tissue samples were preserved in −30 °C at the Center for DNA Taxonomy laboratory, Zoological Survey of India (ZSI), Kolkata, for downstream molecular analysis. The genomic DNA isolation, PCR, purification of the PCR products, and bi-directional sequencing were performed by the standardized protocols (Kundu et al. 2018) Figure 1. The generated sequences were checked through SeqScanner V1.0 (Applied Biosystems Inc., CA, USA), nucleotide BLAST (https://blast.ncbi.nlm.nih.gov/), ORF finder (https://www.ncbi.nlm.nih.gov/orffinder/), and MEGAX (Kumar et al. 2018) to diminish the low-quality sequence read, mismatches, and gaps. The consensus sequences were contributed to the global GenBank database to acquire the accession numbers. To estimate the genetic divergence, phylogenetic analysis, and haplotyping, 13 publicly available sequences of O. hannah and 14 sequences of B. fasciatus were downloaded from GenBank with their collection locality information. Both generated and database sequences were aligned together by using ClustalX to build a dataset (Thompson et al. 1997). The Kimura 2 parameter (K2P) genetic distances for the studied dataset were estimated by MEGAX. The best fit model for the Bayesian analysis (BA) was calculated through Mr.Modeltest v2 with the lowest BIC value (Nylander 2004). The BA phylogeny was constructed in Mr. Bayes 3.1.2 by selecting nst = 6 and rates = invgamma for GTR + G + I model. The MCMC (one cold and three hot chains) was run for 1,000,000 generations with 25% burn-in and trees saving at every 100 generations along with other strand parameters (Ronquist and Huelsenbeck 2003). The constructed BA phylogeny was further beautified in the web-based iTOL tool (https://itol.em) (Letunic and Bork 2007). The database sequence of Lycodon pictus (Accession No. MN395830) under family Colubridae was used as an out-group for the present phylogenetic tree. To understand the genealogical relationships, the haplotype networks were constructed within the different populations of O. hannah and B. fasciatus. The numbers of haplotypes were generated using DnaSP v6 (Rozas et al. 2017). The haplotype diversity (Hd) and nucleotide diversity (p) for the different population of both the species were also calculated through DnaSP v6. The haplotype data files were further used to construct the haplotype network by using PopART (http://popart.otago.ac.nz) (Leigh and Bryant 2015) with the Templeton, Crandall and Sing (TCS) method (Clement et al. 2000).

Figure 1. (A) Live photograph of O. hannah, (B) Live photograph of B. fasciatus, (C) Bayesian phylogeny based on partial mtCytb gene inferred monophyletic clustering of both Elapidae species. (D) TCS network of O. hannah and (E) TCS network of B. fasciatus reveled distinct haplotype of both species collected from Mizoram state in northeast India. Haplotypes are shown in different color circles as represent by collection localities marked in the phylogeny.

Results and discussion

The state of Mizoram considered as a part of Indo-Burma biodiversity hotspot and holds a rich assemblage of herpetofauna (ZSI 2007; Das 2010). The state houses over 60 snake species of which 11 are deadly noxious and cause 2–3 human deaths every year. These venomous snakes are predominately belonging to two families, Elapidae and Viperidae. Among them, the living Kraits are mostly active and bite during night; however, the Vipers and Cobra bite during the day time or dusk. Due to the lack of knowledge on the habitats and characteristics of extant species, snakes are facing man-made threats in this state too. Although the morphological characters readily identify this group of species, their asymmetric body size, coloration, and mimicry frequently led to the doubt in identification (Davis Rabosky et al. 2016). The use of molecular tools not only assist in the rapid species identification but also assure their evolutionary trends, diversification, and other genomic features to discriminate the population (Burbrink and Lawson 2007; Castoe et al. 2012). Due to the efficacy of mitochondrial Cytb gene for determining boundaries between snake species (Slowinski and Keogh 2000; Laopichienpong et al. 2016), the present study aimed to generate the molecular data of two morphologically identified Elapidae species from Mizoram state and contributed to the global databases (B. fasciatus (ZSI_SHT1): GenBank Accession No. MN853157, BOLD Process ID SKBIO001-20 and O. hannah (ZSI_SHT17): GenBank Accession No. MN853158, BOLD Process ID SKBIO002-20). Further, based on the locality information, the database sequences were merged in the analyzed dataset for estimating the intra-species genetic distances, phylogeny and haplotypes.

A total of 50 and 57 variable sites were diagnosed within the studied mitochondrial Cytb dataset (464 bp) of B. fasciatus and O. hannah, respectively. Both the Elapidae species showed 22.1% mean K2P genetic distance with each other. The within-group mean genetic distance of O. hannah was 4.2%, ranging from 0% to 7.6%. All the sequences generated from different localities were clustered together in the present BA phylogeny and a monophyletic clade was formed for this species (Figure 1C) . The studied specimen collected from Mizoram showed close clustering with the specimen (Accession No. EF694840) collected from other parts of India with 2% genetic distance; however, 0.7% genetic distance was revealed with the Burmese specimen (Accession No. AF217842). Further, the northeast Indian specimen showed 2.7% genetic distance with Chinese specimen, 4.1–4.7% genetic distance with Vietnamese specimens, and 7.6% genetic distance with Thailand specimens. The studied sequences of O. hannah resulted in 12 haplotypes with 37 polymorphic sites, haplotype diversity= 0.989. The TCS network showed a distinct haplotype of O. hannah collected from Mizoram (Figure 1D).

The within-group mean genetic distance of B. fasciatus was 3%, ranging from 0% to 8.2%. All the sequences were clustered together in the BA tree and a monophyletic clade was formed. The studied specimen collected from Mizoram showed close clustering with the specimen (Accession No. AF217830) collected from Myanmar with 1.3% genetic distance. Further, the northeast Indian specimen showed 3.3% genetic distance with the specimen collected from Thailand, Vietnam, Thailand, China, and Australia, 3.4% genetic distance with Indonesian specimens, 4% genetic distance with the specimen collected from Indonesia, China, and Thailand. Further, a single database sequence generated from Thailand depicted high genetic distance (6.1%) with the northeast Indian specimen. Based on the low genetic distances of B. fasciatus from different distant geographical localities, the present study assumed that the population of this species might have coalesced within its range distribution. The studied sequences of B. fasciatus resulted in nine haplotypes with 15 polymorphic sites, haplotype diversity = 0.885. The TCS network showed a distinct haplotype of B. fasciatus collected from Mizoram (Figure 1E).

This genetic information is not only adjudicated to the species identification and systematics studies but also triggered other ecological and biological studies of these two venomous snakes from northeast India. The in-depth study of these species will further facilitate the medical emergencies and detect the trade route of snake venom trafficking within and outside of India. We recommend a similar approach can be adopted for other deadly snakes distributed in India to reveal their genetic signature. Apart from this, the awareness among the community, capacity building and training, and improvement of medical facilities in the rural areas could diminish snakebite mortality in northeastern states and other parts of India.

Acknowledgements

We are thankful to the Director of Zoological Survey of India (ZSI), Ministry of Environment, Forest and Climate Change (MoEF&CC), Govt. of India, for providing necessary molecular facilities and support for the study. We are grateful to Michael Vanlalchhuana, Lalrinsanga, Lalmuansanga for their assistance during the field visit. The author (HTL) would like to express his sincere thanks to the Chief Wildlife Warden, Department of Environment, Forest and Climate Change, Govt. of Mizoram, for the collection permission and the DST (SERB), New Delhi for financial support.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

The first author (SK) acknowledges the fellowship grant received from the Council of Scientific and Industrial Research (CSIR) Senior Research Associateship (Scientists’ Pool Scheme) Pool No. 9072-A.