Genetic Diversity and Phylogeography of Thottapalayam thottimvirus (Hantaviridae) in Asian House Shrew (Suncus murinus) in Eurasia

Murid and cricetid rodents were previously believed to be the principal reservoir hosts of hantaviruses. Recently, however, multiple newfound hantaviruses have been discovered in shrews, moles, and bats, suggesting a complex evolutionary history. Little is known about the genetic diversity and geographic distribution of the prototype shrew-borne hantavirus, Thottapalayam thottimvirus (TPMV), carried by the Asian house shrew (Suncus murinus), which is widespread in Asia, Africa, and the Middle East. Comparison of TPMV genomic sequences from two Asian house shrews captured in Myanmar and Pakistan with TPMV strains in GenBank revealed that the Myanmar TPMV strain (H2763) was closely related to the prototype TPMV strain (VRC66412) from India. In the L-segment tree, on the other hand, the Pakistan TPMV strain (PK3629) appeared to be the most divergent, followed by TPMV strains from Nepal, then the Indian-Myanmar strains, and finally TPMV strains from China. The Myanmar strain of TPMV showed sequence similarity of 79.3–96.1% at the nucleotide level, but the deduced amino acid sequences showed a high degree of conservation of more than 94% with TPMV strains from Nepal, India, Pakistan, and China. Cophylogenetic analysis of host cytochrome b and TPMV strains suggested that the Pakistan TPMV strain was mismatched. Phylogenetic trees, based on host cytochrome b and cytochrome c oxidase subunit I genes of mitochondrial DNA, and on host recombination activating gene 1 of nuclear DNA, suggested that the Asian house shrew and Asian highland shrew (Suncus montanus) comprised a species complex. Overall, the geographic-specific clustering of TPMV strains in Asian countries suggested local host-specific adaptation. Additional in-depth studies are warranted to ascertain if TPMV originated in Asian house shrews on the Indian subcontinent.

Murid and cricetid rodents were previously believed to be the principal reservoir hosts of hantaviruses. Recently, however, multiple newfound hantaviruses have been discovered in shrews, moles, and bats, suggesting a complex evolutionary history. Little is known about the genetic diversity and geographic distribution of the prototype shrew-borne hantavirus, Thottapalayam thottimvirus (TPMV), carried by the Asian house shrew (Suncus murinus), which is widespread in Asia, Africa, and the Middle East. Comparison of TPMV genomic sequences from two Asian house shrews captured in Myanmar and Pakistan with TPMV strains in GenBank revealed that the Myanmar TPMV strain (H2763) was closely related to the prototype TPMV strain (VRC66412) from India. In the L-segment tree, on the other hand, the Pakistan TPMV strain (PK3629) appeared to be the most divergent, followed by TPMV strains from Nepal, then the Indian-Myanmar strains, and finally TPMV strains from China. The Myanmar strain of TPMV showed sequence similarity of 79.3-96.1% at the nucleotide level, but the deduced amino acid sequences showed a high degree of conservation of more than 94% with TPMV strains from Nepal, India, Pakistan, and China. Cophylogenetic analysis of host cytochrome b and TPMV strains suggested that the Pakistan TPMV strain was mismatched. Phylogenetic trees, based on host cytochrome b and cytochrome c oxidase subunit I genes of mitochondrial DNA, and on host recombination activating gene 1 of nuclear DNA, suggested that the Asian house shrew and Asian highland shrew (Suncus montanus) comprised a species complex. Overall, the geographic-specific clustering of TPMV strains in Asian countries suggested local host-specific adaptation. Additional in-depth studies are warranted to ascertain if TPMV originated in Asian house shrews on the Indian subcontinent.

Ethics Statement
The guidelines of the American Society of Mammalogists (Kirkland, 1998;Sikes and Animal Care and Use Committee of the American Society of Mammalogists, 2016) were followed for trapping and euthanasia of shrews and for tissue collection and processing. And approvals were obtained from the Ministry of Agriculture and Rural Development in Vietnam and the Institutional Animal Care and Use Committee of the National Institute of Infectious Diseases to conduct the study (permission numbers: 108074, 111126, 112152, 115162, 118180).

RT-PCR and DNA Sequencing
Nested primers for TPMV and other recently identified shrewborne hantaviruses were used to initially screen tissues for hantavirus RNA (Song et al., 2009;Kang et al., 2011c; , 2013). Thereafter, amplification of the full-length S-, M-, and L-genomic segments was attempted. Oligonucleotide primer sequences have been deposited as Supplementary Data S2. Firstand second-round PCR was performed in 20-µL reaction mixtures, containing 250 µM dNTP, 2 mM MgCl 2 , and 0.25 µM of each primer. LA Taq hot start version (Takara Bio) and AmpliTaq gold 360 DNA polymerase (Applied Biosystems, Foster City, CA, USA) were used at 1 U each for the firstand second-round PCR, respectively (Arai et al., 2016b). Initial denaturation at 94 • C for 2 min was followed by two cycles each of denaturation at 94 • C for 30 s, two-degree step-down annealing from 48 to 38 • C for 40 s, and elongation at 68 • C for 1 min, then 32 cycles of denaturation at 94 • C for 40 s, annealing at 42 • C for 40 s, and elongation at 68 • C for 1 min, in a Veriti thermal cycler (Applied Biosystems) and Mastercycler X50 (Eppendorf, Hamburg, Germany) (Arai et al., 2008(Arai et al., , 2012. Amplicons were treated with Exonuclease I and Shrimp Alkaline Phosphatase (New England Biolabs, Ipswich, MA, USA) for 30 min. DNA was sequenced directly using an ABI Prism 377XL Genetic Analyzer (Applied Biosystems) (Arai et al., 2007;Kang et al., 2011c).

Genetic and Phylogenetic Analysis
Partial S-, M-, and L-segment nucleotide and amino acid sequences, amplified from Asian house shrews, were aligned with available hantavirus sequences, using the ClustalW in BioEdit (Thompson et al., 1994). The degree of sequence homology was assessed by pair-wise comparisons (Kang et al., , 2011a. Phylogenetic trees were constructed using MrBayes 3.1.2 (Ronquist and Huelsenbeck, 2003), with the GTR+I+Γ model of evolution, as selected by using jModelTest version 2.1.7 (Darriba et al., 2012). Bayesian analysis consisted of 10 million Markov chain Monte Carlo generations to ensure convergence across two runs of six chains each, with average standard deviations of split frequencies <0.01 and effective sample sizes well over 100, resulting in consensus trees supported by posterior-node probabilities (Kang et al., , 2011a. The co-evolutionary relationships between hantaviruses and their shrew and rodent reservoir hosts were analyzed by the comparative concordance between host and hantavirus cladograms in TreeMap 3b1243 (Charleston and Robertson, 2002;Kang et al., 2009;Arai et al., 2012).

mtDNA and Nuclear Genes Sequencing and Host Phylogeny
To verify the geographic diversity of Asian house shrews and to study their phylogenetic relationships, genomic DNA was extracted from lung tissue using the MagDEA R DNA 200 (GC) (Precision System Science). The entire 1,140-nucleotide cytochrome b (cytb) gene of mitochondrial DNA (mtDNA), the 1,545-nucleotide cytochrome c oxidase subunit I (COI) gene and the recombination activating gene 1 (RAG1) were amplified using the following primer sets: Cy-14724F (5 ′ -GACYARTRRCA TGAAAAAYCAYCGTTGT−3 ′ )/Cy-15909R (5 ′ -CYYCWTYIY TGGTTTA CAAGACYAG−3 ′ ) (Arai et al., 2016b) and KOD multi-enzyme (Toyobo, Osaka, Japan), MammMt-5533F (5 ′ - (Arai et al., 2019) and Phusion enzyme (New England Biolabs), and newly designed primers RAG1-61F (5 ′ -TCTGCACCYGATGAAAT TCARCACC−3 ′ )/RAG1-3139R (5 ′ -CTCCATTGAATCTTGG CTTTCC−3 ′ ) and KOD multi-enzyme, respectively. PCR was performed in 50-µL reaction mixtures, containing 200 µM dNTP and 1 U of KOD multi and Epi DNA polymerase or Phusion enzyme. Initial denaturation was at 95 • C for 2 min, followed by two cycles each of denaturation at 95 • C for 15 s, two-degree step-down annealing from 60 to 50 • C for 30 s, and elongation at 68 • C for 1 min 30 s, then 30 cycles of denaturation at 95 • C for 15 s, annealing at 55 • C for 30 s, and elongation at 68 • C for 1 min 30 s, in a Veriti thermal cycler (Arai et al., 2019). PCR products were purified by Mobispin S-400 (Molecular Biotechnology, Lotzzestrasse, Germany) and were sequenced directly (Arai et al., 2012(Arai et al., , 2019. The models of host nucleotide evolution were selected under jModeltest version 2.1.7, the GTR+I+Γ model for host phylogenetic sequence set, the TrN+G for Cytb, TIM3+ I for COI of and TPM2uf +I models for RAG1 of Suncus sequence sets. The results of modeltest were shown in Supplementary Data 3A-G. Host phylogenetic analysis also consisted of 10 million Markov chain Monte Carlo generations to ensure convergence across two runs of six chains each, with average standard deviations of split frequencies <0.01 and effective sample sizes well over 100, resulting in consensus trees supported by posterior-node probabilities.

Hantavirus Detection
In all but two of the 198 shrew lung tissue samples, multiple attempts to detect hantavirus RNA were unsuccessful ( Table 1). The exceptions were one of 11 and one of three Asian house shrews from Pakistan (captured in Karachi: 24.947802 N,67.122999 E) and Myanmar (captured near a cattle farm in Taung gyi, Shan state: 20.804169 N, 97.060360 E, detail in Supplementary Table 1), respectively, collected in 2013. Sequence analysis of the amplicons revealed TPMV. Amplification of the full-coding region of the S segment and the partial M and L segments was achieved for TPMV strain H2763 (Myanmar), while only partial L-segment sequences were obtained for TPMV strain PK3629 (Pakistan).

Nucleotide and Amino Acid Sequence Analysis
Analysis of the S-, M-, and L-segment sequences of TPMV strain H2763 from Myanmar indicated an overall genomic organization similar to prototype TPMV strain VRC66412 from India. The 1,506-nucleotide S-genomic segment encoded a nucleocapsid (N) protein of 435 amino acids, possibly starting at nucleotide position 68, and a 130-nucleotide 3 ′ -non-coding region. The TPMV S-genomic segment, like that of other recently described hantaviruses in shrews, did not contain the hypothetical NSs open reading frame, typically found in hantaviruses harbored by cricetid rodents.

Hantavirus Phylogeny
TPMV strain H2763 from Myanmar appeared as one cluster in phylogenetic trees, based on the S-, M-, and L-segment sequences, using the Bayesian methods (Figure 2). The TPMV strain PK3629 from Pakistan also constructed one cluster in a tree based on the L segment. The phylogenetic trees suggested that the primordial strain of TPMV originated in northern India and surrounding countries, including Pakistan or Nepal (Figure 2).

Pair-Wise Alignment and Comparison
Pair-wise alignment and comparison of the S segment (1,506 nucleotides), M segment (2,382 nucleotides), and L segment (4,963nucleotides) revealed that TPMV strain H2763 from Myanmar exhibited high sequence similarity to prototype TPMV strain VRC66412 from southern India. The TPMV strain PK3629 from Pakistan showed low nucleotide sequence similarity (79.3%) in the L segment, but the encoded amino acid sequences were highly conserved (94.1-99.2%) with TPMV strains from Nepal, India, Pakistan, and China ( Table 2). Compared with representative hantaviruses from rodents, shrews, and bats, the TPMV strain from Myanmar differed by ∼20-60% at the nucleotide and amino acid levels for each segment.

Co-phylogenetic Analysis of Asian House Shrew and TPMV
As evidenced by co-phylogeny mapping, using a consensus tree based on L-segment sequences, TPMV strains segregated according to the geographic locations of the Asian house shrews (Figures 3B,C). The phylogenetic positions of TPMV strains based on the S and M segments mirrored the phylogenetic relationships of their Asian house shrews, except for the Pakistan strain in the L-segment tree. The Pakistan strain was mismatched between virus and host phylogeography ( Figure 3A).

Phylogenetic Analysis of Asian House Shrew
The molecular identification of TPMV-infected shrews was confirmed as S. murinus murinus by amplification and sequencing of the cytb and COI genes of mtDNA and RAG1 gene of nuclear DNA. Phylogenetic analysis based on the cytb gene indicated Asian house shrews and Etruscan shrews (Suncus etruscus) were clearly distinct (Figure 4). However, the relationships between Asian house shrews, Asian highland shrews (Suncus montanus) and, some S. murinus subspecies (such as S. m. murinus, S. m. kandianus, and S. m. caerulescens) (Meegaskumbura et al., 2010) were less clear. Segregation of Asian house shrews in Asia and Africa was demonstrated by co-phylogeny mapping, using consensus trees based on the cytb, COI and RAG1 genes (Figure 5). The phylogenetic positions of RAG1 in nuclear DNA and cytb and COI in mtDNA were not synchronized for each gene. These data suggest that S. murinus and S. montanus are hybrid species and comprise the S. murinus-S. montanus species complex.

DISCUSSION
The Asian house shrew, one of 18 species in the genus Suncus, is widely distributed throughout Asia and the Pacific, , and M-segments (C), respectively. Letterings for taxa are shown in purple for Sri Lanka, green for Myanmar, blue for Nepal, pink for Bangladesh and India, red for Pakistan, orange for China, black for the other countries, and out groups in both trees. Nepal Asian house shrew was adapted experimental animal strain NP6362. Bangladesh experimental strain BD6364 was used for alternative Indian strain.
Frontiers in Cellular and Infection Microbiology | www.frontiersin.org FIGURE 4 | Bayesian phylogenetic tree, based on the 1,140-nucleotide cytochrome b gene of mtDNA of small mammals within the order Eulipotyphla (families Talpidae and Soricidae), order Rodentia (families Muridae and Cricetidae) and order Chiroptera, suborder Yinpterochiroptera (families Pteropodidae, Hipposideridae, Rhinolophidae), and suborder Yangochiroptera (families Nycteridae, Emballonuridae and Vespertilionidae). The tree was rooted using Elephantulus (order Macroscelidea, GenBank accession numbers DQ901019, DQ901206, and DQ901201) as the outgroup. Numbers at nodes indicate posterior probability values (>0.7) based on 150,000 trees: two replicate Markov chain Monte Carlo runs, consisting of six chains of 10 million generations each sampled every 100 generations with a burn-in of 25,000 (25%). Scale bars indicate nucleotide substitutions per site. Letterings for taxa are shown in green for bats, blue for shrews, purple for moles, black for rodents, red for Elephantulus, and red bold for Asian house shrew. The GenBank accession number for the cytb sequence for Asian house shrews are MT344840 in Nepal and JF784169 in China.
Africa, and the Middle East (Figure 1). It is peridomestic, typically found within areas of human habitation, and has become dependent on discarded human food waste. Asian house shrews may have been intentionally introduced by humans, similar to Rattus rodents, into Africa (Egypt, Eritrea, Kenya, Republic of Djibouti, Rwanda, Sudan, and Tanzania), the Middle East (Iran, Iraq, Kingdom of Bahrain, Kuwait, Saudi Arabia, Sultanate of Oman, and Yemen), the islands within the Indian Ocean (Comoros, Republic of Madagascar, Republic of Mauritius, and Réunion), and Asia and the Pacific (Japan, Guam, and Philippines) (Kang et al., 2011c). Genetic analysis and treemap dendrograms of RAG1 and COI, RAG1 and cytb, and COI and cytb suggest that Asian house shrews may represent hybrids with the Asian highland shrew in Sri Lanka and some area of Eurasia (Figure 5).
The previously held conventional view that hantaviruses coevolved with their reservoir hosts has been challenged recently by the conjecture that preferential host switching and local hostspecific adaptation account for the congruent phylogenies of hantaviruses and their small mammal hosts (Ramsden et al., 2009). Multiple examples of host sharing are now known for hantaviruses hosted by rodents and shrews . Whether or not TPMV exhibits such host sharing with evidence of carriage by other species of the genus Suncus requires future investigation.
Based on phylogenetic analysis of mtDNA and nuclear genes, as well as karyotype and morphological analysis, the taxonomy FIGURE 5 | Comparisons of nuclear and mitochondria genes generated by TreeMap 3, using Bayesian method, based on the recombination activating gene 1 (RAG1), cytochrome b (cytb) gene, and cytochrome oxidase I (COI) gene of the Asian house shrew and related shrews. (A) Left tree on the RAG1 sequence, while the mitochondria tree on the right was based on the COI gene or (B) cytb gene, and (C) also comparison of cytb and COI genes, respectively. The nuclear (RAG-1) and mitochondria relationship (cytb and COI) were listed in Supplementary Table 1. Letterings for taxa are shown in right table. The deposited cytb sequences are attached dot. The tree was rooted using Suncus etruscus (order Eulipotyphla, family Soricidae, GenBank accession numbers: cytb, LC126597; COI, MK410384; and RAG1, MT344767). The RAG1 order in Figures 5A,B was the same. The order of COI and cytb in Figure 5C was adjusted for optimal.
of the Asian house shrew is still unclear. Asian house shrews include at least two subspecies (S. murinus murinus, S. murinus kandianus, and S. murinus caerulescens), and the Asian highland shrew (S. montanus) is morphologically very similar. Our genetic analysis suggests that morphological based S. m. murinus, S. m. kadianus, S. m. caerulescens, and S. montanus represent hybrid species. Thus, a species complex has been proposed (Ohdachi et al., 2016).
Although genetically diverse strains of TPMV have been detected in Asian house shrews from Nepal (Kang et al., 2011c) and China (Guo et al., 2011), the geographic distribution and evolutionary origins of TPMV are still unclear. Our data FIGURE 6 | Map of distribution of Asian house shrew and phylogenetic tree based on L segment of TPMV. (A) TPMV cluster based on L-segment in phylogenetic analysis. (B) Blue arrows in map were estimated expansion root based on TPMV phylogeny. Pakistan strain was captured at red star, Nepal strains (Kang et al., 2011c) were captured at light blue circle, Indian strain (Carey et al., 1971) was captured at pink circle, Myanmar strain was captured at green star and Chinese strains (Guo et al., 2011) were captured at red circles. Pakistan and Myanmar strains were collected in this study (star symbols). Pink area is shown Asian house shrew distribution (Ohdachi et al., 2016). suggest the possibility that TPMV expanded from the Indian subcontinent (Figure 6). The evolutionary time scale of TPMV is faster than that of its host and the host is older than the ancient trade routes between the Middle East and China.

CONCLUSION
Disappointingly, TPMV was detected in lung tissues of only two Asian house shrews, one from Myanmar and one from Pakistan. The reasons for this are not entirely clear, but it might be the result of the focal nature of TPMV infection, as is typical of other hantaviruses. Future studies on the phylogeography of TPMV and the Asian house shrew should provide valuable insights into the geographic radiation.

ETHICS STATEMENT
The animal study was reviewed and approved by National Institute of Infectious Diseases (NIID), Institutional Animal Care and Use Committee.