Evolution of Hemoglobin Genes in a Subterranean Rodent Species (Lasiopodomys mandarinus)

The Mandarin vole (Lasiopodomys mandarinus), a typical subterranean rodent, has undergone hematological adaptations to tolerate the hypoxic/hypercapnic underground environment. Hemoglobin (Hb) genes encode respiratory proteins functioning principally in oxygen binding and transport to various tissues and organs. To investigate the evolution of α- and β-hemoglobin (Hb) in subterranean rodent species, we sequenced Hb genes of the Mandarin vole and the related aboveground Brandt’s vole (L. brandtii). Sequencing showed that in both voles, α-globin was encoded by a cluster of five functional genes in the following linkage order: HBZ, HBA-T1, HBQ-T1, HBA-T2, and HBQ-T2; among these, HBQ-T2 is a pseudogene in both voles. The β-globin gene cluster in both voles also included five functional genes in the following linkage order: HBE, HBE/HBG, HBG, HBB-T1, and HBB-T2. Phylogenetic analysis revealed that the Mandarin vole underwent convergent evolution with its related aboveground species (Brandt’s vole) but not with other subterranean rodent species. Selection pressure analyses revealed that α- and β-globin genes are under strong purifying selection (ω < 1), and branch-site analyses identified positive selection sites on HBAT-T1 and HBB-T1 in different subterranean rodent species. This suggests that the adaptive evolution of these genes enhanced the ability of Hb to store and transport oxygen in subterranean rodent species. Our findings highlight the critical roles of Hb genes in the evolution of hypoxia tolerance in subterranean rodent species.


Introduction
Subterranean rodent species spend their whole lives in complete darkness at relatively stable temperatures in hypoxic and hypercapnic underground burrow systems [1,2]. Although hypoxic and hypercapnic environments present major challenges for the survival of subterranean rodent species, previous reports have noted that Spalax carmeli, which is indigenous to Israel, inhabits heavy clay soil and can survive under minimal oxygen (O 2 ) levels (7.2%) and maximal carbon dioxide (CO 2 ) levels (6.1%) [3] and that the naked mole-rat (Heterocephalus glaber), which is a burrowing rodent native to parts of East Africa, can tolerate hours of extreme hypoxia and survive for 18 min under total O 2 deprivation (anoxia) without any apparent injury [4]. The ability of subterranean rodent species to survive under conditions of high CO 2 and low O 2 levels without severe deleterious physiological effects or behavioral changes [4][5][6] suggests that it will be interesting to investigate the specific adaptations of mammals that live in naturally rodent species.
Mechanisms underlying hypoxia tolerance identified in different subterranean rodent species include modified blood properties, anatomical and biochemical changes in respiratory organs, and differences in the structures and functions of gene products [7][8][9]. The important adaptation of subterranean rodent species to hypoxic environments are the increased blood vessel density, increased erythrocytosis, and improved hemoglobin (Hb) capacity of O 2 delivery and storage [6,8,10]. Hb is responsible for the transport of O 2 from the lungs to O 2 -demanding tissues, suggesting that this protein plays an important role in the survival of subterranean rodent species under low ambient O 2 conditions. Hb is a heterotetramer comprising two α-chain (HbA) and two β-chain (HbB) subunits that form two semirigid dimers (α1β1 and α2β2) [11,12]. The functional adaptability of Hb genes manifest both convergent and divergent adaptive evolutionary characteristics via recombination events, unequal gene exchange, genetic conversion, and mutations that aid in long-term survival in chronically hypoxic environments [13][14][15]. For example, intraspecific allele substitution involving multiple polygenic combinations has been described in the deer mouse (Peromyscus maniculatus) [16]; β15Trp has been shown to change to Phe in the fetal domestic yak (Bos grunniens) and β135Ala to Val in the adult yak [14]. Further, a positive selection site has been reported at position144 of the HbB subunit of gymnotiform electric fish [17], and parallel β13 and β83 variants have been documented in Andean hummingbirds [18]. All of these genetic alterations increase the affinity of Hb for O 2 . Although studies have investigated hypoxia-induced physiological adaptation in subterranean rodent species with a focus on the roles of vascular endothelial growth factor (VEGF), hypoxia-inducible factor-1α (Hif-1α), erythropoietin (EPO), myoglobin, neuroglobin, and cytoglobin [8], there is limited research on the adaptive evolution of Hb in subterranean rodent species.
The Mandarin vole (Lasiopodomys mandarinus) is a type of subterranean rodent species capable of surviving in environments with minimal O 2 levels (16.04%) and maximal CO 2 levels (2.55%) [19]. In the present study, we cloned and sequenced Hb genes from the Mandarin vole and the related aboveground Brandt's vole (L. brandtii). The main objectives of the present study were to assess the structure and function of Hb in subterranean rodent species adapted to hypoxia and to provide indirect evidence of the physiological function of Hb in subterranean rodent species to highlight the biomedical significance of this protein.

Hb Gene Structure
Sequencing showed that the α-globin gene cluster in the Mandarin vole included four functional globin genes: one embryonic globin gene (HBZ), two adult globin genes (HBA-T1 and HBA-T2), and one θ globin gene (HBQ-T1). In the Mandarin vole, HBQ-T2 is a pseudogene, and in terms of chromosomal structure beginning with HBZ and terminating with HBQ-T2. The β-globin gene in the Mandarin vole is a cluster of five functional globin genes in the following linkage order: one embryonic globin gene (HBE), one chimeric globin gene (HBE/HBG), one globin gene (HBG), and two adult globin genes (HBB-T1 and HBB-T2) (Figure 1). The αand β-globin genes in Brandt's vole have the same structure as those in the Mandarin vole, with their amino acid compositions in the two tested species being slightly different, coding sequences and differences in the number of amino acids encoded by Hb genes in the two vole species are shown in Table 1.

Phylogenetic Trees
The topology of Maximum likelihood (ML) trees was similar to that of Bayesian trees, and Bayesian trees were selected to present the final results. Bayesian trees of αand β-globin genes were both constructed using HKY+G as the ideal model. The topology of their Bayesian trees indicated a lack of convergent evolution in subterranean rodent species, the Mandarin vole, Brandt's vole, and the prairie vole (Microtus ochrogaster), all of which belong to the Microtinae subfamily within the Cricetidae Biology 2020, 9, 106 3 of 10 family, were clustered together. The Damaraland mole-rat (Fukomys damarensis) and naked mole-rat, which belong to the Bathyergidae family, were clustered together; however, Spalax was located in a separate single branch (Figures 2 and 3).

Phylogenetic Trees
The topology of Maximum likelihood (ML) trees was similar to that of Bayesian trees, and Bayesian trees were selected to present the final results. Bayesian trees of α-and β-globin genes were both constructed using HKY+G as the ideal model. The topology of their Bayesian trees indicated a lack of convergent evolution in subterranean rodent species, the Mandarin vole, Brandt's vole, and the prairie vole (Microtus ochrogaster), all of which belong to the Microtinae subfamily within the Cricetidae family, were clustered together. The Damaraland mole-rat (Fukomys damarensis) and naked mole-rat, which belong to the Bathyergidae family, were clustered together; however, Spalax was located in a separate single branch (Figures 2 and 3).

Selection Pressure Analysis
Branch models were used to calculate the ω-value of α-and β-globin genes; the estimated ω ratio of the α-and β-globin genes was <1, indicating that most Hb genes had undergone purifying selection (Supplementary Table S1). To investigate whether there are restricted positive selection sites in some specific lineages, positive selection sites were identified using the branch-site model for Hb genes ( Table 2)

Discussion
The performance of the hematological system of subterranean rodent species is well established. Hb plays important roles in increased O 2 transport in response to hypoxia [12,16]. Hb gene structure and adaptive evolution have been extensively studied in species living in hypoxic environments at high altitudes and in deep ocean [20,21]. However, the Hb gene structure and evolution in subterranean rodent species remains comparatively unclear. Hence, the objective of the present study was to elucidate the evolutionary adaptation of Hb genes in subterranean rodent species to hypoxic environments.
The organization of the αand β-globin gene clusters in the Mandarin vole parallels that in Brandt's vole, i.e., the structure and number of functional copies of the β-globin gene clusters in both voles are similar, as documented in the deer mouse and prairie vole [22]. The structure and number of functional copies of the α-globin gene cluster were consistent between the two tested voles but differed in other rodent species. Therefore, Hb polymorphism exists among different rodent species, but no significant differences were observed in this regard between the two tested voles despite their different living environments.
Phylogenetic reconstruction of Hb genes revealed that the Mandarin vole did not undergo convergent evolution with other subterranean Spalax, the naked mole-rat, or the Damaraland mole-rat, but it underwent convergent evolution with the related aboveground Brandt's vole and prairie vole. The latter two species and the Mandarin vole belong to the Microtinae subfamily, and their phylogeny has been shown to follow that of traditional rodent taxonomic groups [23,24]; therefore, we hypothesized that Hb genes in the Mandarin vole underwent conservative evolution.
Determination of ω ratios revealed that all Hb genes have undergone purifying selection, indicating that the Hb protein has an important function in O 2 transport. However, the results of our analyses support the hypothesis that HBA-T1 and HBB-T1 in different subterranean rodent species have been subjected to positive selection, suggesting the evolutionary optimization of functionality under different environmental conditions. For example, in HBA-T1, amino acid substitutions in S. galili were located at a proximal binding site situated near or within the glycosylation and metal-binding domains responsible for cotranslational and post-translational protein modifications and O 2 binding, respectively [21]. In naked mole-rat, substitution of residues 87 and 89 was close to the proximal histidine, which is implicated in the iron-proximal histidine linkage and seems to be essential in the maintenance of Hb-O 2 binding [17]. Further, the substitution of residue 124 in the H helical region may change the spacing of hydrogen bonds [25], thereby affecting the molecular structure of proteins.
Furthermore, positive selection sites are located in the A, E, G, and H helical regions of HBB-T1 in S. galili, the naked mole-rat, and the Damaraland mole-rat, and they may be essential for homology-based structural models [12]. Site 143 in the H helical region of HBB-T1 in the Damaraland mole-rat is adjacent to β-147 histidine. A previous study showed that His-HC3 (147) of β-Hb in the common carp (Cyprinus carpio) plays a key role in the root effect, which is associated with non-cooperative O 2 binding and decreased O 2 affinity [17]. Another study also demonstrated that β-globin underwent adaptive evolution in the A helical region in the pika; therefore, selected 17 sites in S. galili in their study may have important Hb-associated functions in O 2 binding [25]. Other positive selection sites were localized to residues of postulated regions and may affect these functions [26][27][28]. Thus, we hypothesized that positive selection sites in these two genes are involved in modulating Hb-O 2 binding and we demonstrated adaptability to different hypoxia-associated elevational zones in two subterranean species.
It is worth mentioning that there were no noteworthy sites on Hb genes in the vole mandarin and Brandt's vole. This is possibly because there is no O 2 deficiency in the living environment of Brandt's vole as it is an aboveground rodent. Although the Mandarin vole is a subterranean rodent, the O 2 level in their burrow system is relatively low compared with that in the environments of other subterranean species [3,4,19]. For example, S. galili can survive under conditions of minimal O 2 levels (12.7%) and the naked mole-rat can tolerate hours of extreme hypoxia and survive 18 min of total O 2 deprivation (anoxia) without any apparent injury. The lowest O 2 level measured in high-water-content soil of a burrow system of the Mandarin vole was only 16.04%, but the physiological characteristics of blood in this species including increased HIF-1α and VEGF levels and blood vessel density and decreased EPO level under hypoxia, indicating adaption to hypoxia.

Animals
Juvenile Mandarin voles were captured on a farmland located in Xinzheng, Henan Province, China (34 •

DNA Extraction, Amplification, and Sequencing
Total DNA was extracted from the brain tissues using the TIANamp Genomic DNA Kit [TIANGEN Biotech (Beijing) Co., Ltd., Beijing, China]. DNA quality was analyzed by agarose gel electrophoresis; DNA concentration was assessed using the NanoDrop ND-1000 spectrophotometer (NanoDrop Technologies, LLC, Wilmington, DE, USA). DNA used for amplification of Hb genes had integral bands, and the A260/A280 or A260/A230 ratios ranged from 1.9 to 2.1.
Hb genes were amplified by polymerase chain reaction (PCR) using primers designed with Oligo 7.0 (https://www.oligo.net$\delimiter"026E30F$) with reference to sequences from closely related species (Table 4). PCR products were separated by electrophoresis using 0.1% agarose gels. Single-purpose fragment DNA strips were cut out from agarose gels, purified, and recovered using the QIAquick Gel Extraction Kit (Beijing ComWin Biotech Co., Ltd., Beijing, China). Complete sequences were inserted into the pMD19-T vector (Takara Bio, Inc., Shiga, Japan). After bacterial selection and colony PCR verification, the Sanger method was performed to sequence the PCR products using GENEWIZ (Suzhou, China). After combining the genomic information of the Mandarin vole, the complete Hb gene structure was diagramed.

Data Analysis
The Hb gene coding nucleotide sequences of other species were downloaded from the GenBank database (https://www.ncbi.nlm.nih.gov/genbank/) and aligned using the Basic Local Alignment Search Tool (http://www.ncbi.nlm.nih.gov). The sequence of each gene was aligned using the MUSCLE algorithm. α-globin/β-globin gene clusters were combined using SequenceMatrix 1.7.8 for the construction of phylogenetic trees [29]. Their Bayesian trees were constructed using MrBayes 3.2.6 [30,31], and the best-fit model for the MrBayes tree was selected using Modeltest. Finally, four Markov chains were run for 2 × 10 6 generations and sampled every 1000 generations for Bayesian tree construction. ML trees for the best nucleotide models were constructed using MEGA × 10.1, and confidence levels of tree topologies were assessed using 1000 bootstrap replicates.
We determined selection pressure for each Hb subunit gene using CODEML, a part of PAML 4.7 [32]. First, we used branch models to calculate the ω-value (dN/dS = nonsynonymous/synonymous) of Hb genes, with ω < 1, ω = 1, and ω > 1 indicating purifying, neutral, and positive selection, respectively. We used the one-ratio model that enforces the same ω ratio for all lineages and compared it with ω = 1 (i.e., ω = 1 for all lineages) to check the selection pressures of the two genes.
Second, we used a branch-site model (Model A, A null) to analyze positive selection sites on specific branches of each gene. Then, in PAML, we performed LRT to determine the differences in log-likelihood values between the two nested models using χ2 distribution, with an aim to determine whether the differences were significant by comparing the nested models. Finally, site models M7 and M8 were used to compute recombination sites in all assessed species [33], and the GARD method within Data Monkey (https://www.datamonkey.org/) was performed to identify the number and location of breakpoints and sequences involved in putative recombination events in different genes [34].

Conclusions
Hb genes are evolutionarily conserved and play important functions in O 2 storage and transport in subterranean rodent species. Elucidation of the evolution of Hb genes offers a valuable framework for the identification and study of the mechanism underlying adaptation to hypoxic environments in subterranean rodent species. Further functional experiments including functional characterization of Hb proteins are expected to completely reveal these mechanisms.
Author Contributions: Conceptualization, Z.W.; investigation, K.Y., D.L., D.P. and S.G.; analysis, writing-original draft preparation, and writing-review and editing, H.S. All authors have read and agreed to the published version of the manuscript.