Continued Co-circulation and Evolution of the Ancient Subcontinent Lineage Despite Predominance of the Arctic-Like Lineage Rabies Viruses (RABV) in India

Keywords

Rabies viruses • Evolution • Arctic like lineage • Subcontinent lineage

Introduction

Rabies is one of the longest known human diseases and was first documented at least 4,000 years ago [1-3]. Rabies is among the key Neglected Tropical Diseases (NTDs) that predominantly affects poor and vulnerable populations. World Health Organization (WHO) in 2010 estimated that rabies was responsible for 26,400 to 61,000 human deaths annually and a majority of these deaths (84%) occur in rural areas [4].

Rabies is caused by Rabies Virus (RABV), which belongs to the genus Lyssavirus of the family Rhabdoviridae and order Mononegavirales [5]. RABV virions are enveloped, rod-or bullet-shaped, with single stranded negative sense RNA genomes measuring approximately 12 kb. Among the 12 identified species within the genus Lyssavirus, RABV is the most widely distributed geographically and infect a wide range of animals within the orders Carnivora and Chiroptera [6-8]. Rabies is an acute and highly fatal disease characterized by encephalomyelitis. All terrestrial mammals, including humans are known to be susceptible to rabies. Saliva of infected animals is an important source of RABV and is transmitted to susceptible animals either by bite, scratch or mucous membrane exposure. Although there are effective post-exposure vaccines that can prevent rabies, rabies still causes approximately 55,000 human deaths annually with a majority (>95%) of these occur in Asia and Africa [9]. According to 2011 WHO estimates, 1.74 million disability-adjusted life years are lost each year due to rabies with an estimated annual economic burden of US$ 583.5 million. While rabies can be prevented by vaccination, currently there is no treatment after the onset of the clinical disease.

In most developing nations, domestic dog continues to be the key reservoir of rabies and plays an important role in interspecies transmission of RABVs. In addition, several wild carnivorous species are known to maintain RABVs in limited geographic regions with a very limited contribution to the overall human rabies burden. While many livestock species such as cattle are susceptible to rabies, they only act as epidemiological dead-end hosts and do not transmit the disease further [10]. Multiple lineages of dog RABVs have been described that include African, Asian, Arctic-like, Cosmopolitan and Indian subcontinent.

Rabies remains a major public health problem in India and a WHO sponsored national rabies survey in 2006 found that 20,000 human deaths due to rabies occur each year in India [11]. In 2011 molecular epidemiology based on nucleoprotein (N) and Phosphoprotein (P) genes found that all the Indian RABV isolates belonged to the Arctic/ Arctic-like viruses and clustered within a single clade [12]. Subsequently, a study in 2014 using the ecto-domain coding region of the glycoprotein gene sequence, found co-existence of both Arctic like 1 lineage and Indian subcontinent lineage RABV in India [13]. However, a study in 2015 reported that based on G gene sequence, all the Indian isolates obtained between 2001-2014 from six different species of animals were genetically related to Arctic-like 1a lineage viruses [14]. Indian RABV isolates clustered within the Arctic/Arcticlike clade have been found to significantly diverge evolutionarily from the cosmopolitan and other lineages of RABVs circulating across Southeast Asia [15]. Most of the RABV molecular epidemiological studies carried out in India were done with isolates collected over a long period of time and from multiple geographical locations across the country.

Molecular epidemiologic approaches are very helpful to track the spread of variant RABV into new geographic regions [16]. Further, molecular epidemiological analyses of RABVs in Rabies endemic countries such as India would allow accurate analysis of the spread and evolution of RABVs to design better control methods. A previous study found a RABV isolate from the city of Chennai (previously known as Madras) in the southeastern coast in India clustered more closely to a variant RABV found in Sri Lanka. Chennai, in the state of Tamil Nadu, is the Southernmost city in India. The movement of humans and their animals between Sri Lanka and India was proposed to be a source of the emergence of the variant RABV in India. Movement of people continues to occur between the southeastern coastal area of India and Sri Lanka. To better understand the genetic divergence of RABV lineages circulating in India and to investigate the proposed RABV incursions from Sri Lanka, we carried out molecular evolutionary analyses of RABV isolates collected in 2014 from Southern Indian states of Tamil Nadu and Goa.

Materials and Methods

Samples and rabies testing

Twenty-six brain samples from animals died of rabies were used in this study. The samples included twenty-three isolates from domestic dogs, two samples from goats and one sample from cattle. The samples used in this study were collected in 2014 from two different states on the East and West coast of southern India (Tamil Nadu and Goa) (Table 1). The samples were tested at the rabies unit of Madras veterinary college for rabies confirmation by Fluorescent Antibody Test (FAT) [17]. In addition to FAT, N-gene specific RT-PCR [18] using 20%brain homogenate was carried out to detect rabies viral genomes. Sample collection from animals was performed according to the guidelines approved by the Institutional Animal Ethics Committee (IAEC) of Tamil Nadu Veterinary and Animal Sciences University (TANUVAS), Chennai, India.

Table 1. List of rabies positive samples and their origin.

RNA extraction and RT-PCR

Extraction of total RNA from the brain tissue homogenates was carried out using TRIzol reagent (Invitrogen), following the manufacturer’s instructions. N and G gene specific primers were designed using the available rabies sequences of the Indian isolates on public databases (Table 2).

Table 2. Gene specific primers used for gene amplification and sequencing.

G gene was amplified by RT-PCRs as two overlapping fragments. Using 1 μg of total RNA complementary DNA (cDNA) synthesis was done using superscript reverse transcriptase III enzyme (Invitrogen) per manufacturer’s instruction. PCR amplification was done using platinum Pfx polymerase (Invitrogen) and thermal cycler profile was as follows: polymerase activation at 94°C for 5 min; 35 cycles of denaturation at 94°C for 15 s, annealing at 55°C for 30’s and extension at 68°C for 1 min; extension at 72°C for 10 min.

Sequencing

Amplified gene segments were purified by QIAquick gel purification kit (Qiagen) per the manufacturer’s instructions. DNA sequencing of the N and G genes were performed using the Illumina HiSeq2000 platform (Illumina). For 26 isolates of RABV, complete nucleotide sequencing of Nucleoprotein (N) gene (1353 nucleotides) and partial sequencing of glycoprotein (G) gene (711 nucleotides) was performed.

Phylogenetic and evolutionary analysis

N and G-gene sequences belonging to global RABV isolates available as of Dec 2019 were retrieved from GenBank and used for phylogenetic analysis. Multiple sequence alignment was carried out using MAFFT [19,20]. Recombination detection was carried out using RDP4. Phylogenetic reconstruction was carried out using both alignment-based and alignment-free methods. Alignment-based phylogenetic tree was derived using the Maximum Likelihood (ML) method available in the PhyML package with 1000 bootstrap replicates.

Molecular clock analysis was carried out using BEAST v1.10.4. Earlier model selection was carried out using jModelTest and bat rabies virus was used as an outgroup. Molecular clock behavior was tested using TempEst. Relaxed clock with lognormal distribution was used for molecular clock analysis with GTR+I+gamma as substitution model and constant coalescence as demographic model. Markov Chain Monte Carlo (MCMC) algorithm was run for 10 billion steps and sampled every 10000 steps. Tracer v1.7 was used for assessing convergence and iTOL server was used for visualization of the phylogenetic trees.

Results and Discussion

Phylogenetic and evolutionary analysis based on N gene

Sequencing of nucleoprotein gene of 26 Indian isolates was carried out, of which 20 sequences were of complete length (1353 nucleotides). MVC-29 (GenBank: MH258816, Chennai, dog), MVC-31 (GenBank: MH258818, Chennai; dog), MVC-2 (GenBank: MH258835, Goa; dog), MVC-3 (GenBank: MH258834, Goa; dog), MVC-5 (GenBank: MH258832, Goa; dog) and MVC-45 (GenBank: MH258822, Chennai; goat) are observed to share 100% identity.

The 20 nucleoprotein sequences of Indian isolates along with global data of RABV sequences (totaling 389 entries; Supplementary file 1) were used for Molecular Phylogenetic Analysis (MPA). Multiple sequence alignment revealed that sequences shared ~31% identity and ~45% similarity. The maximum likelihood phylogenetic tree showed clustering based on known lineages, namely Arctic-like, Africa 2, Indian subcontinent, Cosmopolitan and Asian. Molecular-clock behavior was absent when data pertaining to bat-RABV was included in the analysis. Hence, the mean rate of nucleotide substitution of N gene belonging to dog-RABV was estimated to be 2.69 × 10^-4 substitutions/site/yr (95% Highest Posterior Density (HPD) 2.19 × 10^-4 to 3.34 × 10^-4). Based on the N gene, the RABV isolates sequenced in this study clustered into two lineages, Arctic-like (15 isolates) and Indian subcontinent (5 isolates). Indian RABV isolates belonging to Arctic-like lineage reported in this study clustered with sub-lineage Arctic-like lineage 1a (AL-1a) which contains few members from neighboring countries like Nepal and Afghanistan. The Indian subcontinent lineage RABV isolates reported in this study clustered with previously reported RABV isolates from India and Sri Lanka belonging to Indian subcontinent lineage.

Phylogenetic tree of nucleoprotein (N) sequences of 15 Rabies Virus (RABV) Indian isolates along with RABV sequences from across the globe (totaling 389 sequences) derived using Maximum likelihood method available in PhyML with 1000 bootstrap replicates. The Indian isolates sequenced in this study that belong to Arctic-like and Indian sub-continent lineages have been highlighted in green and yellow respectively.

Phylogenetic and evolutionary analysis based on G gene

Sequencing of RABV glycoprotein gene was carried out for 26 Indian isolates, of which 22 sequences of comparable length (711 nucleotides) were analyzed. Similar to the observation in analysis of complete N gene, partial G-gene sequences of isolates MVC-4 (GenBank: MH258811, dog), MVC-29 (GenBank: MH258800, dog), MVC-19 (GenBank: MH258807, dog), MVC-30 (GenBank: MH258799, dog) and MVC-38 (GenBank: MH258796) are exactly identical. MVC-2 (GenBank: MH258813, dog), MVC-3 (GenBank: MH258812, dog) and MVC-15 (GenBank: MH258808, dog) are identical. MVC-21 (GenBank: MH258815, dog), MVC-34 (GenBank: MH258802, dog) and MVC-23 (GenBank: MH258803, dog) are identical. MVC-49 (GenBank: MH258793, cattle), MVC-48 (GenBank: MH258806, dog) and MVC-51 (GenBank: MH258805, dog) are identical.

The 22 glycoprotein sequences of Indian isolates along with global data of RABV sequences (total 374 entries; Supplementary file 1), were used for molecular phylogenetic analysis. MSA showed ~38% identity and ~65% similarity. As observed in N-gene analysis, molecular clock behavior was absent upon inclusion of bat RABV sequences. Hence, the mean rate of nucleotide substitution estimated from the partial glycoprotein sequences for dog-RABV was 4.5 × 10^-4 substitutions/site/yr (95% HPD 3.44 × 10^-4 to 5.56 × 10^-4). Of the 22 Indian RABV isolates sequenced in this study, 17 RABV clustered with Arctic-like lineage (sub-lineage AL-1a) and 5 clustered with Indian subcontinent lineage.

The phylogenetic trees based on the N and G genes displayed similar topologies, indicating the presence of equivalent clades in both trees. For example, both trees showed Indian subcontinent lineage RABV isolates clustering closely with other Indian subcontinent lineage RABV isolates from Sri Lanka. Notably, ML tree showed that the Arctic like RABV isolates from this study clustered closely with human RABV isolates from India and Germany.

Phylogenetic tree of partial glycoprotein (G) sequences of 21 Rabies virus (RABV) Indian isolates along with RABV sequences from across the globe (totaling 374 sequences) derived using Maximum likelihood method available in PhyML with 1000 bootstrap replicates. The Indian isolates sequenced in this study that belong to Arctic-like and Indian sub-continent lineages have been highlighted in green and yellow respectively.

The time-scaled evolutionary trees obtained using Bayesian analysis (data not shown) based on both N and G genes indicated that the isolates of Arctic-like lineage diversified more recently as compared to isolates of Indian subcontinent lineage.

Conclusion

A key aspect of RNA virus evolution is the exceptionally high rate of nucleotide substitution, that is responsible of for the highly variable rates of molecular evolution among viral species. Like all RNA viruses, RABVs constantly evolve and the rate of evolutionary change in different host species is determined by the nature of virus-host interactions. We investigated the evolutionary rates of RABV isolates from India using partial sequences of the Glycoprotein (G), the surface protein that allows RABV to enter the cells of nervous system and the complete gene sequences of Nucleoprotein (N), which is the constituent of viral capsid and involved in virus replication.

The mean rate of nucleotide substitution estimated from the N gene sequences was 2.69 × 10^-4 substitutions/site/yr (95% Highest Posterior Density (HPD) 2.19 × 10^-4 to 3.34 × 10^-4). Interestingly, these estimates are in close agreement with the reported N-gene mean rate of substitution in RABVs from bats and terrestrial mammals suggesting that the dog RABVs in India continue to evolve at the same rate as bat RABVs from around the world. The mean rate of nucleotide substitution estimated from the partial G gene sequences is 4.5 × 10^-4 substitutions/site/yr (95% HPD 3.44 × 10^-4 to 5.56 × 10^-4) which is marginally higher than that reported in a previous study, which has evaluated substitution rate in Indian RABV isolates based on the ecto domain of glycoprotein. Earlier study using the complete G gene sequences reported the mean rate of substitution to be 3.9 × 10⁻⁴ substitutions/site/year (95% HPD 1.2-6.5 × 10⁻⁴). The results of this study suggest a higher mean rate of evolutionary change in G gene in Indian dog RABV isolates as compared to evolutionary rates of lyssaviruses. This difference could be due to use of partial versus full G gene sequences for analysis, hence further in-depth whole genome analyses are warranted to further validate this observation.

The phylogenetic trees constructed using both N and G genes showed RABV isolates from this study belonged to two distinct linages namely Arctic lineage and Indian subcontinent lineage. This finding is in agreement with an earlier study that showed coexistence of two distinct lineages of RABV in India. An earlier study reported genetic clustering of RABV isolates based on the geographical region within India. However, the present study found that RABV isolates clustered based on membership to respective lineages rather than geographic proximity. Notably, this study reports that the RABV isolates belonging to Indian subcontinent lineage clustered closely with the isolates from Sri Lanka that belong to other subcontinent lineage. This highlights the continued incursion and/or circulation of the variant subcontinent lineages of RABVs in India that might have been originated from Sri Lanka. Phylogenetic analysis also revealed that the Arctic like RABV isolates from this study clustered closely with human RABV isolates from India and Germany and they belonged to AL1a sub-lineage described earlier by Troupin, et al.

Complex mechanisms are involved in the maintenance of rabies viruses in the primary host species rather than that involved in serial transmission to a new host species. RABV causes lethal infection in dogs with an infectious period of less than 1 week (and typically only 2–4 days). Hence, RABV need to be constantly transmitted among members of the same species to be maintained at the population. Continued circulation of the two distinct linages of RABVs in India despite several vaccination and stray dog population control campaigns suggest that there is a perpetuated transmission among dogs and enzootic viral establishment of the two distinct lineages among the dog population.

Several studies investigating RABV biology highlighted that RABVs are sensitive to control measures. The maintenance of RABV is suggested to be influenced by a delicate balance between densitydependent transmissions and virus-induced mortality. However, other factors including demographic and spatial structure have recently been suggested to generate observed epidemic cycles for RABV maintenance in a country.

Most of the rabies deaths in India involve children who are more likely to encounter infected stray dogs. Government of India has implement a “National Rabies Control Program” approved during 12^th five-year plan (2012-2017) which aimed to prevent the human rabies deaths and to prevent transmission of rabies through canine (dog) rabies control. Dog Population Management (DPM) involves improving the health and well-being of stray or community dogs by vaccination and reduce dog population size by routine birth control programs which can facilitate more effective rabies control. Despite vaccination for many years, there was no downward trend in the RABV incidence in India and there is continued co-circulation of two distinct lineages as found by this study. Intentional or unintentional translocation of dogs between different geographical regions could potentially compromise natural or vaccine-generated barriers. As the Indian subcontinent lineage RABV isolates from Chennai are closely related to viruses from Sri Lanka, this suggests the possibility of regular incursions of RABVs between the two countries. Hence, it is also important to check the cross-boundary moment of dogs and implement proper quarantine to prevent the spread of rabies.

Acknowledgments

The authors would express their sincere gratitude to Late Dr. Elankumaran Subbiah who played a key role in the conception and design of these experiments, who sadly could not be a co-author of this publication. This study was funded by a UKIERI Trilateral Research in Partnership Award (ND/CONT/E/12-13/704). The authors would also like to acknowledge the access to computational facilities at the bioinformatics centre, SPPU, supported by the department of biotechnology, Govt. of India.

References

1. GM, BAER. "Rabies-an historical perspective." Inect Agent Dis 3 (1994): 168-180.
2. Belotto A, LF Leanes, MC Schneider and H Tamayo, et al. "Overview of rabies in the Americas." Virus Res 111 (2005): 5-12.

   [Crossref] [Google Scholar] [PubMed]

3. Fisher, Christine R, Daniel G Streicker and Matthias J. Schnell. "The spread and evolution of rabies virus: conquering new frontiers." Nat Rev Microbiol 16 (2018): 241-255.

   [Crossref] [Google Scholar] [PubMed]

4. Dacheux, Laurent, Olivier Delmas and Hervé Bourhy. "Human rabies encephalitis prevention and treatment: progress since Pasteur's discovery." Infect Disord Drug Targets 11 (2011): 251-299.

   [Crossref] [Google Scholar] [PubMed]

5. Nel, Louis H and Wanda Markotter. "Lyssaviruses." Crit Rev Microbiol 33 (2007): 301-324.

   [Crossref] [Google Scholar] [PubMed]

6. Nel LH and CE Rupprecht. "Emergence of lyssaviruses in the Old World: the case of Africa." Curr Top Microbiol Immunol 315 (2007): 161-193.

   [Crossref] [Google Scholar] [PubMed]

7. Knobel, Darryn L, Sarah Cleaveland, Paul G. Coleman and Eric M. Fèvre, et al. "Re-evaluating the burden of rabies in Africa and Asia." Bull World Health Organ 83 (2005): 360-368. [Crossref] [Google Scholar]

   [PubMed]

8. Nel, Louis H, Louise H Taylor, Deepashree Balaram and Kim AS Doyle. "Global partnerships are critical to advance the control of neglected zoonotic diseases: the case of the global alliance for rabies control." Acta tropica 165 (2017): 274-279.

   [Crossref] [Google Scholar] [PubMed]

9. Sudarshan MK, SN Madhusudana, BJ Mahendra and NSN Rao, et al. "Assessing the burden of human rabies in India: results of a national multi-center epidemiological survey." Int J Infect Dis 11 (2007): 29-35.

   [Crossref] [Google Scholar] [PubMed]

10. Reddy, GB Manjunatha, R Singh, RP Singh and KP Singh, et al. "Molecular characterization of Indian rabies virus isolates by partial sequencing of nucleoprotein (N) and phosphoprotein (P) genes." Virus Genes 43 (2011): 13-17.

    [Crossref] [Google Scholar] [PubMed]

11. Reddy, RV Chandrasekhar, B. Mohana Subramanian, KSNL Surendra and RP Aravindh Babu, et al. "Rabies virus isolates of India–simultaneous existence of two distinct evolutionary lineages." Infect Genet Evol 27 (2014): 163-172.

    [Crossref] [Google Scholar] [PubMed]

12. Cherian, Susan, Rajendra Singh, KP Singh and GB Manjunatha Reddy, et al. "Phylogenetic analysis of Indian rabies virus isolates targeting the complete glycoprotein gene." Infect Genet Evol 36 (2015): 333-338.

    [Crossref] [Google Scholar] [PubMed]

13. Nadin-Davis, Susan A, Geoff Turner, Joel PV Paul and Shampur N. Madhusudana, et al. "Emergence of Arctic-like rabies lineage in India." Emerg Infect Dis 13 (2007): 111.

    [Crossref] [Google Scholar] [PubMed]

14. Nadin-Davis, Susan A, G. Allen Casey and Alexander I. Wandeler. "A molecular epidemiological study of rabies virus in central Ontario and western Quebec." J Gen Virol 75 (1994): 2575-2583.

    [Crossref] [Google Scholar] [PubMed]

15. Aravindh Babu, RPS Manoharan, P Ramadass and NDJ Chandran. "Evaluation of RT-PCR assay for routine laboratory diagnosis of rabies in post mortem brain samples from different species of animals." Indian J Virol 23 (2012): 392-396.

    [Crossref] [Google Scholar] [PubMed]

16. Edgar, Robert C. "MUSCLE: a multiple sequence alignment method with reduced time and space complexity." BMC Bioinformatics 5 (2004): 1-19.

    [Crossref] [Google Scholar] [PubMed]

17. Edgar, Robert C. "MUSCLE: multiple sequence alignment with high accuracy and high throughput." Nucleic Acids Res 32 (2004): 1792-1797.

    [Crossref] [Google Scholar] [PubMed]

18. Martin, Darren P, Carolyn Williamson and David Posada. "RDP2: recombination detection and analysis from sequence alignments." Bioinformatics 21 (2005): 260-262.

    [Crossref] [Google Scholar] [PubMed]

19. Guindon, Stéphane and Olivier Gascuel. "A simple, fast and accurate algorithm to estimate large phylogenies by maximum likelihood." Syst Biol 52 (2003): 696-704.

    [Crossref] [Google Scholar] [PubMed]

20. Drummond, Alexei J, Marc A. Suchard, Dong Xie and Andrew Rambaut. "Bayesian phylogenetics with BEAUti and the BEAST 1.7." Mol Biol Evol 29 (2012): 1969-1973.

    [Crossref] [Google Scholar] [PubMed]