- Research article
- Open Access
Adaptive patterns in the p53 protein sequence of the hypoxia- and cancer-tolerant blind mole rat Spalax
BMC Evolutionary Biology volume 16, Article number: 177 (2016)
The subterranean blind mole rat, Spalax (genus Nannospalax) endures extreme hypoxic conditions and fluctuations in oxygen levels that threaten DNA integrity. Nevertheless, Spalax is long-lived, does not develop spontaneous cancer, and exhibits an outstanding resistance to carcinogenesis in vivo, as well as anti-cancer capabilities in vitro. We hypothesized that adaptations to similar extreme environmental conditions involve common mechanisms for overcoming stress-induced DNA damage. Therefore, we aimed to identify shared features among species that are adapted to hypoxic stress in the sequence of the tumor-suppressor protein p53, a master regulator of the DNA-damage response (DDR).
We found that the sequences of p53 transactivation subdomain 2 (TAD2) and tetramerization and regulatory domains (TD and RD) are more similar among hypoxia-tolerant species than expected from phylogeny. Specific positions in these domains composed patterns that are more frequent in hypoxia-tolerant species and have proven to be good predictors of species’ classification into stress-related categories. Some of these positions, which are known to be involved in the interactions between p53 and critical DDR proteins, were identified as positively selected. By 3D modeling of p53 interactions with the coactivator p300 and the DNA repair protein RPA70, we demonstrated that, compared to humans, these substitutions potentially reduce the binding of these proteins to Spalax p53.
We conclude that extreme hypoxic conditions may have led to convergent evolutionary adaptations of the DDR via TAD2 and TD/RD domains of p53.
The blind mole rat of the genus Nannospalax (hereafter, Spalax) is a solitary subterranean mammal  that experiences extreme and abrupt fluctuations in O2/CO2levels. Spalax survives low-oxygen content (~7 % O2) in its natural underground habitat, and even lower (3 % O2) under laboratory conditions [2–4]. One of the greatest challenges faced by Spalax is repeated exposure to acute hypoxia followed by rapid re-oxygenation, which leads to oxidative stress . Hypoxia and oxidative stress are two types of cellular stressors that place a great risk on cellular functions and genomic stability [6–8]. Importantly, the combination of these stressors seems to further enhance genomic instability compared to the separate effects of each stress type. Hypoxia depletes dNTPs and represses DNA repair pathways that are required to overcome replication-stress, while re-oxygenation induces replication-restart that coincides with reactive oxygen species (ROS)- induced oxidative DNA damage, at a time when the repair pathways have not yet recovered [9, 10].
Genomic instability underlies both cancer and aging . Nevertheless, Spalax does not develop spontaneous tumors and does not show clear age-related phenotypic changes, despite its relatively long lifespan (~20 years in captivity) . Moreover, Spalax displays an outstanding tolerance to chemically induced carcinogenesis in vivo, and its fibroblasts inhibit cancer growth in vitro . Another long-lived (~30 years) subterranean rodent in which spontaneous cancer is an uncommon phenomenon, relative to similar-sized rodents, is the naked mole rat, Heterocephalus glaber (hereafter, H. glaber) [14, 15]. H. glaber is phylogenetically distant from Spalax (estimated divergence time is 77.9 million years) . Yet, similar to Spalax, it shows adaptations to hypoxic stress  and in vitro ability to inhibit cancer cell growth . Bearing in mind the well-known strong positive correlation between aging and cancer epidemiology in many metazoans, including humans , common features in Spalax and H. glaber deserve special attention. Such features may relate to evolutionary adaptations to common environmental conditions that might play a compensatory role in dealing with consecutive stress-induced DNA damage.
The subterranean habitat is not the only environment in which oxygen fluctuations may occur. Diving mammals experience rapid transitions from apnea to re-oxygenation , while hibernating mammals face reduced metabolism and fluctuations in blood flow and oxygen consumption [20–23]. “Metabolic shutdowns” happen also in desert mammals, such as the desert mouse  and the jerboa . Evidence for oxidative stress or antioxidant adaptations were indeed found in Spalax , diving mammals [19, 27], hibernating mammals , and desert mammals . Bats are another example of mammals showing antioxidant adaptations , as might be expected from the high metabolic stress these mammals experience during flight [30, 31]. Interestingly, bats’ lifespan increases with hibernation, body mass, and the occasional use of caves . Their exceptional longevity is explained by multiple mechanisms for resisting oxidative damage  and is attributed to positively selected genes in DNA-damage-checkpoint pathways . The above mentioned studies indicate that environmental conditions, which include extreme changes in oxygen supply, may require adaptations to cope with oxidative stress that are not restricted to antioxidant activity, but also include changes in the DNA-damage response (DDR) as a second line of defense.
Sophisticated DDR pathways constantly monitor genome integrity  and delay or stop cell-cycle progression at critical stages in response to unrepaired DNA damage, thereby preventing replication of the damaged DNA . The tumor suppressor protein p53 interacts with multiple proteins participating in complex pathways of the DDR . It plays an important role in cellular pathways that control genomic instability and is involved in critical defense and regulatory pathways such as cell-cycle arrest, senescence, and apoptosis . The tight control by senescence and apoptosis prevents the uncontrolled proliferation of damaged cells; however, it may also deplete stem and progenitor cell pools, thus promoting tissue degeneration and aging . This emphasizes the subtle role of p53 in regulating pathways related to cancer and aging . Spalax has acquired an R174K amino acid substitution in the DNA-binding domain (DBD) of p53. Mutation at this position in human p53 leads to impaired induction of apoptotic and cell-cycle arrest genes and is known as a frequent mutation found in many human cancers [39, 40]. The fact that Spalax does not develop spontaneous cancer implies that there are changes in other regions of the p53 protein and/or other related proteins to ensure integrity of the related signaling pathways and to maintain homeostasis.
Due to the critical role of p53 in the regulation of apoptosis and DNA repair, we aimed to investigate whether Spalax p53 protein sequence shares structural features with other phylogenetically distant mammals that are also adapted to acute-fluctuating hypoxia. Such features may indicate joint adaptations of the cellular stress response via p53 that are related to hypoxic stress resistance. Here, we have identified changes in p53 domains, which are known to include intrinsically disordered regions  that are highly diverse in evolution , and found that these domains' sequences are similar among mammals adapted to stressful hypoxic conditions, more than expected from phylogeny. These domains harbor binding sites of p53 with proteins that participate in the DDR and in certain metabolic pathways such as replication protein A 70-kDa DNA-binding subunit (RPA70), the histone acetyltransferase p300, and silent information regulator 2 (Sir2) [43–45]. Thus, changes in these domains could modulate the way in which p53 orchestrates these pathways under extreme hypoxic conditions in the subterranean habitat.
Different intra-group similarities for hypoxia-tolerant and hypoxia-sensitive taxa
We hypothesized that certain functional domains in p53 have convergently evolved among hypoxia-tolerant species. If this is true, it is expected, for example, that some domains in Spalax p53 are more similar to the corresponding domains of other hypoxia-tolerant (HT) species, such as diving mammals, than to those of hypoxia-sensitive (HS) species, such as murine rodents, although Spalax is much closer to murines. To test this, we used a balanced experimental design (Fig. 1a) in which a single representative HT species and a single representative HS species is tested for different taxonomic groups (e.g., carnivores, rodents, even-toed ungulates, insectivores, etc.). In public databases we were able to find five representative HT species whose p53 genes were sequenced, belonging to five different mammalian taxonomic orders. Accordingly, five pairs of closely related species were used, where in each pair one is HT and the other is HS. In this test, amino acid sequence distances were calculated according to the Poisson correction model (see Methods) for the following p53 domains (Fig. 1b): transactivation subdomain 1 (TAD1), transactivation subdomain 2 (TAD2), proline-rich domain (PRD), DNA binding domain (DBD), nuclear localization signal (NLS), and tetramerization and regulatory domains (TD/RD). For each domain in p53, it was then tested whether the pairwise distances within the group involving HT species are smaller than those of the matched group of HS species (H1 hypothesis). In addition, the test was repeated using the estimated pairwise divergence time published in the literature .
Following the test described in Methods, we calculated Mann–Whitney statistics U that characterizes the difference in intra-group similarities (based on all-to-all distances in each group). P-values for this statistics were calculated in two ways: using U distribution (Mann-Whitney ) and distribution obtained by Monte-Carlo simulations (Table 1). Monte-Carlo simulations were used to avoid biases caused by intra-group dependencies derived from the usage of each species more than once when calculating the distances of all possible pairs of species in each group (see Methods). Monte-Carlo simulations were conducted with k HT = k HS = 5 and n HT = n HS = 5*4/2 = 10, where k represents the number of species in each group (five in our test) and n represents all possible pairs of species in each group (ten in our test). U-statistic for TAD2, DBD, and TD/RD domains was significantly smaller (with FDR adjusted p < 0.05), allowing the rejection of the H0 hypothesis that components of vectors S HT and S HS have the same joint distribution (see Methods). These observed significant differences are unlikely to be explained only by phylogeny, as no significant differences in the intra-group domain similarities were detected using the estimated pairwise divergence time.
To further verify that the obtained results remain consistent when a wider range of species is used (see Methods), we also compared Spalax p53 protein sequence, as a reference sequence to: (1) the entire NCBI nr database; (2) p53 sequences of 46 other species, which include all HT species whose p53 sequence was sequenced, as well as additional species that are adapted to other types of metabolic stress (Additional file 1: Table S1, Additional file 2: Figure S1, hypoxic- metabolic- and non-stress groups). The results of these tests (Additional file 3: Table S2) were consistent with the results based on the balanced experimental design for TAD2 and TD/RD (but not DBD), and indicated that these domains in Spalax are more similar to those of HT species than to those belonging to species with the highest phylogenetic relatedness to Spalax.
Stress-related patterns in p53 TAD2, TD, and RD
The results above may indicate that unique sequence patterns have convergently evolved in some p53 domains of hypoxia-tolerant species. Accordingly, we investigated candidate sequence patterns in p53 that could specify the similarity between species adapted to stressful hypoxic conditions. Here, “pattern” means the same combination of residues, in the same order, with a constant number of residues between them. We have identified two types of patterns in TAD2 and TD/RD domains of p53: (1) patterns that are more frequent in species from the hypoxic-/metabolic-stress groups relative to species from the non-stress group; (2) patterns that are characteristic solely to all “hypoxic-stress” species (Additional file 1: Table S1) and exist in none of the species included in the metabolic-/non-stress groups.
To identify the first type of patterns, the complete p53 sequences of 47 species from the three above-mentioned groups were aligned by MAFFT (sorted by input order that was determined according to the respective groups). Positions that led to incorrect matching between homologous residues were manually corrected. Common patterns in TAD2 and TD/RD were manually identified and the frequencies of these patterns in the hypoxic-/metabolic-stress groups relative to the non-stress group were subjected to Fisher’s 2 × 2 exact test. In addition, the significance of these patterns in all three groups was tested by Fisher’s 2 × 3 exact test. This data set was used as it represents a wide range of species in terms of stress adaptations, so that (1) we could verify that the patterns identified in p53 protein sequence exist in most or all of hypoxia-tolerant species (2) patterns identified as different between the two extremes (hypoxic vs. non-stress groups), could be further evaluated in mammals experiencing other types of stresses (metabolic stress group). Nevertheless, as this data set is not absolutely phylogeneticaly balanced (see Methods), only the identified patterns for which the p-value was smaller than 0.005 were selected and retested in the mentioned balanced experimental design (see Fig. 1a).
Among the identified patterns (Fig. 2) that were found to be significantly more frequent in species included in the hypoxic-stress group relative to the non-stress group (Additional file 4: Table S3), some of the patterns were also found in the metabolic-stress group, but with a lower level of significance. In TAD2, the pattern LLXXE was significantly more frequent in the hypoxic (p < 0.0001) and metabolic (p < 0.05) groups compared to the non-stress group. Position A52 (Spalax p53; corresponds to E51 in human p53), and its related patterns EXXA and LLXXEXXA, were significantly more frequent only in the hypoxic-stress group (p < 0.005 and 0.0001, respectively). In TD/RD, the patterns KXXGE and LM were significantly more frequent in the hypoxic-stress group (p < 0.0001 for both patterns) and in the metabolic-stress group (p < 0.005 and 0.05, respectively) relative to the non-stress group. The patterns LLXXEXXA, EXXA and LM in TAD2 and RD (but not in TD) were still significantly more frequent in the hypoxic-stress group vs. the non-stress group of the balanced experimental design (p < 0.05), despite the small number of species included in each group (Additional file 4: Table S3, Additional file 5: Figure S2).
The next step of our analysis was to check whether the residues participating in the patterns identified in TAD2, TD, and RD (see shaded columns in Fig. 2a) are sufficiently informative to enable the robust classification of species into the three groups: hypoxic-stress, metabolic-stress, and non-stress. We used the “classification trees” algorithm to classify an extended list of mammals (66 rather than the initial 47 species, see Additional file 1: Table S1, Additional file 6: Figure S3), based on the identified positions in p53 as predictors (see Methods). The classification (Additional file 4: Table S3, Additional file 7: Figure S4) showed that all species from the hypoxic-stress group (subterranean and diving mammals) had very high probabilities of being assigned to this group. An analogous result (namely, a classification that corresponds to the subjective assignment) was obtained for the species that belonged to the non-stress group. In the metabolic-stress group, most of the species showed high probability of being assigned to stress-related groups (mostly the metabolic-stress group). Interestingly, in this category, mammals with large body mass tended to have lower probabilities of being classified into stress groups. Nineteen (19) additional species (that are not included in any of the hypoxic/metabolic/non-stress groups, see Methods) were approximately equally distributed between non-stress and stress-related groups. The most important predictor (Additional file 4: Table S3) was found to be p53 position 51 (human p53) in TAD2. These results provide further support that adaptation to stressful hypoxic conditions is predicated by changes at the identified positions.
A second type of pattern in TAD2 (Fig. 2b) domain defines unambiguously (one-to-one) the species included in the hypoxic-stress group. This pattern can be described as a combination of two motifs: [DE]LLX[ST] and [DE]XX[AV]XWL. The combination of these two motifs was found in all hypoxia-tolerant species, but was not found in any of the other species included in the sample (total 47 species). We then checked whether this pattern exists in the 19 additional species (out of the total 66, see Additional file 1: Table S1). The pattern was not found in any of these species, with one exception: the African elephant. In other words, among 66 species tested, this pattern was found only in hypoxia-tolerant species and in the African elephant. Nevertheless, for the African elephant, the pattern in p53-RD (LM) was not found, and thus it is distinguished from the hypoxic-stress group.
In addition, as TAD2 of p53 was shown to resemble TAD2 domain of Krueppel-like factor 1 (KLF1) , we also compared the sequences of KLF1-TAD2 domain of hypoxia-tolerant species to those of hypoxia-sensitive species (according to the balanced experimental design). We did not identify any hypoxia-specific pattern (Additional file 5: Figure S2) in KLF1-TAD2. This may indicate that the considered patterns are associated with unique functionality of p53 rather than with a functionality related to a global hypoxic response.
3D models of p53 in complex with other DDR proteins provide evidence of a unique DDR in spalax
To investigate the possible functional role of the residues identified in p53 as ‘stress-related’ (Fig. 2), we first allocated the available structures of p53 in complex with other proteins in which the interaction includes these residues. We subjected TAD2 and TD/RD sequences to protein BLAST, querying the Protein Data Bank (PDB) database. Only a few structures of the proteins that interact with p53 via these domains have been resolved for human proteins to date (for main examples see Additional file 8: Table S4). We then examined the literature related to these structures, and found that the residues included in the stress-related patterns are indeed directly involved in the interactions of p53 with some of the proteins (Table 2). For example, the pattern ‘EXXA’ (Spalax p53) corresponds to residues ‘DXXE’ (positions 48–51, human p53), which are responsible for electrostatic interactions between p53 and other proteins such as RPA70 and p62.
In addition, we observed large differences between Spalax and humans with respect to p53 interactions with proteins involved in DDR pathways via TAD2, as demonstrated in Additional file 9: Figure S5. Therefore, we hypothesized that p53’s interactions with these proteins are different in Spalax compared to humans. We then specifically examined p53’s interaction with RPA70 by model replacement of the human p53 residues participating in the direct interface with RPA70, with the corresponding residues in Spalax. From this analysis, it seems that the human p53–RPA70 interaction interface involves more hydrophobic and ionic interactions relative to Spalax and therefore, the interaction in Spalax is less stable (Fig. 3a). Importantly, the residues that mediate the interaction with p53 are identical in humans and Spalax in the RPA70 protein (data not shown), unlike the interface-residues in the p53 protein, which are different between the two species. This suggests that p53 and RPA70 did not co-evolve and strengthens the conclusion that the interaction between them is different in Spalax compared to humans. Similar results were obtained when we modeled Spalax-p53’s interaction with p300 via TAD2. The interface in the p300 side of the interaction was conserved (data not shown); however, we have identified reduced electrostatic interactions in Spalax compared to humans due to different interface-residues in p53 (Fig. 3b). We then used MutaBind server to estimate the putative changes in the binding affinity of Spalax p53 interactions with other DDR proteins (Additional file 8: Table S4). This analysis supported the possibility that Spalax’s stress-related substitutions in TAD2 reduce the binding affinity of p53 to other DDR proteins such as RPA70 and p62, compared to humans. For the stress-related pattern in TD/RD (owned by both Spalax and humans) it was demonstrated that this pattern potentially increases the binding affinity of p53 to proteins such as sir2, compared to other species that do not own the pattern.
Site-specific positive selection in p53
The adaptive role of the sites in which residues composing the above-mentioned stress-related patterns are located was further investigated. We hypothesized that specific p53 amino acid sites in TAD2 and TD/RD were positively selected during evolution. By ‘positive selection’, we mean amino acid substitutions that conferred adaptive advantages under selective pressure. PAML package and Selecton online-server were used to estimate codons’ non-synonymous vs. synonymous substitution rates (dN/dS). A likelihood ratio test (LRT) was conducted comparing two nested models: a null model that assumes no positive selection and an alternative model that assumes positive selection. The same 65 sequences were also analyzed using the HyPhy model MEME (Mixed Effects Model of Episodic Diversifying Selection), which can identify sites that were positively selected only in a small fraction of the tree branches. Positions E49, A52 in Spalax p53-TAD2 (corresponding to D48, E51 in the human p53) and position L381 in Spalax p53-RD (corresponding to L383 in the human p53), which participate in the stress-related patterns, were significantly identified as positively selected by 4, 3, and 1 models out of 4, respectively (see Fig. 4a, Additional file 10: Table S5). The putative role of these positions is illustrated in the 3D structures of p53’s interactions with several proteins via TAD2 and RD (Fig. 4b). Accordingly, since TAD2 seems to have a unique stress-related structure and function, positions E49 and A52 (Spalax p53) inside this domain are promising candidates for further analysis. Notably, all tests implemented in Selecton, PAML, and Hyphy identified another residue in the DBD of p53 as a significantly positively selected site: P127 in Spalax p53, corresponding to A129 in human p53. This position is located in a loop that fine-tunes the angle between two stapled beta-strands (Fig. 4b). The importance of this position could be related to the role of the two beta-strands’ edges. For example, these edges include positions such as K120 (human p53), which is directly involved in p53 interactions with the DNA and the anti-apoptotic protein Bcl-xL .
In this study, we identified similarities in p53 functional domains that are shared between mammals adapted to stressful hypoxic conditions. TAD2 and TD/RD of p53 are rapidly evolving, compared to other domains in p53, such as the DBD . Nevertheless, they are similar between mammals that are not phylogenetically closely related, but share the ability to endure extreme changes in oxygen levels. We suggest that this similarity is an indication of convergent evolution, resulting from selective pressures that induced common molecular adaptations via these p53 domains. We further identified sequence patterns in these domains that are more frequent in mammals adapted to acute-fluctuating hypoxia and demonstrated that the residues participating in the identified patterns can predict species classification according to stress characteristics. An additional pattern in TAD2 was observed in all hypoxia tolerant-species and in none of the other species participating in this study (a total of 66 tested), with one exception: the African elephant. Interestingly, this mammal is cancer-resistant and expresses a unique DNA damage response that is attributed to mechanisms related to p53 . The latter pattern is a combination of two motifs: [DE]LLX[ST], which includes Ser46 and may relate to the binding of proteins that induces its phosphorylation, and [DE]XX[AV]XWL, which includes a motif that is parallel to a motif found in p53-TAD1, which is attributed to the binding of p53 with MDM2  and thus may relate to the competitive binding of MDM2 to p53-TAD2 and p53-TAD1 .
Protein sequence motifs, which are located in rapidly evolving intrinsically disordered regions of a hub protein and include segments that mediate its binding with other proteins, fits well with the definition of Short Linear Motifs (SLiMs) [53–55]. These motifs represent a functional change in a flexible low-affinity binding region, which may divert the regulation of a complex signaling network. We suggest that the patterns identified in this study may potentially serve as SLiMs that have evolved among hypoxia-tolerant species. The evolution of SLiMs seems to be very complex and is not yet well understood. Nevertheless, recent theories based on evidence from the past decade  suggest that the appearance of SLiMs could be random and “evolutionarily transient” (i.e., ex nihilo lost and gained along lineages). This unique evolutionary process may minimize the dependencies between phylogenetically related species and between species and their ancestors; however, it does not contradict events of positive and negative selection along their evolution. In this study, we also showed that positions of residues that make up the identified patterns (E49 and A52 in Spalax p53, corresponding to D48 and E51 in human p53), can be considered as positively selected. Importantly, these residues stabilize electrostatic interactions between p53 and proteins such as the cofactor and histone acetyltransferase p300 and the DNA-repair proteins RPA70 and Tfb1/p62 [57–61].
The intrinsically disordered nature of TAD2 and RD enables them to interact with a range of proteins that function as a complex network. Among them are DNA-repair proteins (Additional file 8: Table S4) that when defected cause susceptibility to a variety of congenital cancers and progeroid syndromes (genetic disorders that mimic physiological aging). These syndromes are attributed to the failure of genome-maintenance mechanisms , suggesting the involvement of TAD2 and TD/RD in p53-mediated regulation of genomic stability . Only a few 3-D structures of the proteins in complex with p53 via TAD2 and TD/RD have been solved to date, while other interactions are still unknown. Thus, the possible role of residues making up the patterns identified in this study has yet to be elucidated. Nevertheless, to get indications of the possible adaptive importance of the specific changes in these domains in Spalax p53, we applied these changes in silico in human p53 and showed, albeit indirectly, some evidence for decreased affinity between the modified p53 (according to the hypoxia-tolerant pattern) and RPA70 compared to human p53. This complies with a previous study  demonstrating that a few combinations of p53 mutations: [L43A, M44A, L45A], or [D48H, D49H] or [I50A, E51A, Q52A], have disrupted the binding of p53 to RPA70. Each of these combinations includes residues that participate in the pattern identified in the current study as stress-related in TAD2 (highlighted in bold). Other studies [64, 65] have demonstrated that the mutations D48H and D49H in p53-TAD2 preserved the transactivation function of p53, but disrupted p53 interaction with RPA70 and the subsequent inhibition of homologous recombination (HR) by p53 . This finding suggests the involvement of these residues in the sequestration of RPA70, which is required for HR suppression, and therefore it may be concluded that an inhibition of HR suppression (i.e., HR promotion) might occur in Spalax.
TAD2 of p53 functions as a single-stranded (ss) DNA mimic and interacts with proteins that also interact with ssDNA [58, 59]. This DNA mimicking ability may have a role in the dynamics of these proteins’ translocation to the damage site. For example, RPA70 has higher affinity to ssDNA than to p53, and upon DNA damage, it detaches from p53 and integrates with ssDNA at the damage site . Thus, the potential decrease in the binding affinity of p53 to RPA70 in Spalax may lead to faster recruitment of RPA70 in the case of DNA damage. Another protein that interacts with p53 via TAD2 and shows similar structural competition with ssDNA is breast cancer 2 protein (BRCA2), which is also involved in double-strand break repair and HR [66, 67]. There is evidence that residues in TAD2, including the above-mentioned D48 and E51 (human p53), play a role in the binding of p53 to BRCA2 , implying a similar role in the recruitment of BRCA2 to the damage site. In this interaction, it was suggested that charged residues, such as D48 and E51 (corresponding to E49 and A52 in Spalax p53) attract the partner domains of BRCA2 to p53 at the initial phase of binding, whereas hydrophobic residues, such as L43 and L45 (human p53), which is also included in the identified patterns, are associated with the folding of the participating activation domain in BRCA2.
P53 transactivation is regulated by p300, which also affects p53 turnover, depending on the cellular context and the environmental stimuli, e.g., those that induce DNA damage . A recent study investigating the structure of p53-TAD2 in complex with p300  indicates the involvement of almost all the residues identified in the TAD2-stress-related pattern in this interaction. As in the case of RPA70, mutations in positions D48 and E51 (human p53) were shown to induce a decrease in the binding affinity of p53 to p300 . This fits our modeling of the p53-p300 interaction via TAD2 indicating a potentially weaker binding in Spalax relative to humans. Importantly, under hypoxia, p53 crosstalks with the hypoxic response transcription factor hypoxia inducible factor 1 (HIF-1) and both compete for p300 as a co-activator to enhance transcription of their target genes . This competition may suggest that the possible decrease in the binding affinity of Spalax p53 to p300 indirectly increases p300 availability to HIF-1 and thus enhances HIF-1-mediated transcription (hypoxic response) on the account of p53-mediated-transcription (apoptotic response). This could serve Spalax as a complementary mechanism to its constitutively high level of HIF-1 expression , as the upregulation of HIF-1 alone is not sufficient to enhance its target-genes transcription, and a compatible amount of available p300 is needed. The potential enhancement of HIF-1-mediated transcription on the account of p53 transcription may shift the balance in Spalax cells between survival and apoptosis, in favor of survival.
Most interesting is to examine proteins that interact with both TAD2 and TD/RD of p53, as such interactions can point to a functional link between these two domains. A functional link between TAD2 and RD of p53 was previously demonstrated in the case of the homeodomain-interacting protein kinase-2 (HIPK2), which its knockdown inhibited S46 phosphorylation in TAD2 and K382 acetylation in RD; both S46 phosphorylation and K382 acetylation are required for p53–mediated apoptosis . CREB-binding protein (CBP) belong to the same protein family  as p300 and share a similar structure and function as co-activators. Both interact with p53 via TAD2 and acetylate the same residue (K382, human p53-RD) [73, 74]. This residue is also a substrate for the NAD-dependent protein deacetylase SIRT1 (human sirtuin) . Hence, CBP/p300 and SIRT1 may compete for binding p53 in this location in RD. It was shown that the adjacent L383 (human p53), which is included in the stress-related patterns identified in this study, participates in p53 interaction with CBP and SIRT1 via RD, which is required for the acetylation/deacetylation of K382 in p53 (Table 2, Fig. 4). Acetylation of p53 by p300 positively regulates p53 activity, while deacetylation by SIRT1 negatively regulates p53 activity, suppresses apoptosis, and prolongs cellular survival in response to DNA damage . It can be suggested that the patterns we found in p53, may promote p53 deacetylation by SIRT1 on the account of p53 acetylation by p300/CBP, due to the potential decrease in p300/CBP binding to p53 (via TAD2) that may indirectly increase the availability of K382 in RD to SIRT1. SIRT1 homologs have been shown to induce slow-aging in Caenorhabditis elegans, Drosophila, and mice . It was suggested that this effect is achieved by delaying apoptosis via p53 and giving cells more time to repair the damage . We see the similarities in p53 protein sequence found between Spalax and other stress-tolerant and phylogenetically distant mammals as a gateway for investigating potential evolutionary-convergent interconnected mechanisms that may relate its tolerance to extreme changes in oxygen levels to longevity and cancer resistance.
In this study we showed that the p53 domains TAD2 and TD/RD are more similar between species adapted to hypoxic stress than expected from phylogeny. Patterns shared in these domains among hypoxia-tolerant species include residues that mediate the binding of p53 with critical DDR proteins. Part of the residues are located in positions that were identified as positively selected and potentially play a role in reducing the binding of p53 with p300 and RPA70, suggesting a possible change in p53-mediated regulation of the DDR among species adapted to extreme changes in oxygen levels. Further investigation of TAD2 and TD/RD in the context of hypoxia and DNA repair may lead to progress in understanding cellular hypoxia tolerance, cancer resistance, and aging.
Data set and experimental design
In the current study, p53 sequences of 66 placental mammals were used (Additional file 1: Table S1, Additional file 6: Figure S3). From this group of 66 species, 47 species (Additional file 2: Figure S1) were classified into 3 groups: (1) “hypoxic-stress” group, which includes 15 species adapted to acute and transient environmental hypoxia (e.g., subterranean and diving mammals). Species in this group acquired the ability to inhabit environments that are characterized by extreme and frequent oxygen fluctuations and oxygen deprivation. In this group, we included all known-to-date mammals that, to the best of our knowledge, fit this criterion and have the full length of the p53 protein sequence published. The information regarding mammalian adaptations to different types of hypoxic stresses is available in a review by Ramirez et al. , and in our previous study on Spalax ; (2) “non-stress” group, which includes 16 species that do not endure any type of hypoxic stress (organismal/external/chronic/acute) or other types of stresses affecting cellular metabolism and extreme changes in availability of oxygen to tissues and cells; (3) “metabolic stress” group, which includes 16 species that fit neither the acute-hypoxic nor the non-stress criteria, and are known to endure other types of metabolic stresses such as hibernation, dehydration, extreme temperatures, chronic hypoxia, flying, or different combinations of these stresses (Additional file 1: Table S1). Extra caution should be taken when interpreting the results of the statistical analysis related to this group due to its heterogeneity. The remaining other 19 mammals could not be classified into any of the groups due to insufficient or inconclusive information in the literature or because including them in one of the groups strongly disturb the balance between groups. These mammals were termed as “additional” group and their p53 sequences were mostly used in analyses for which prior classification is not required, such as positive selection.
Species classification into the three first groups specified above was based on the summarized literature regarding each species from the viewpoint of stress adaptations. The classification was optimized as explained below, so that stress adaptations related to p53 could be investigated. For example, some mammals were classified as “non-stress” despite the fact they have semi-aquatic family members (e.g., fishing cat and mink are relatives of cat and ferret, respectively) or evidence of semi-aquatic adaptations (e.g., rhinoceros and humans) . In addition, in this experimental design, the ability to obtain balanced proportions of species from different taxonomic orders and families in each group was also limited by the number of species for which the p53 sequence had been published and by the actual adaptations of the species. Thus, for example, in the metabolic-stress group, there are more rodents than in the other groups. Due to the mentioned limitations, an additional experimental design was used in which each hypoxia-tolerant species was paired with its closest hypoxia-sensitive phylogenetic relative that has an available p53 sequence (see Fig. 1a). Specifically, we selected all pairs of species so that in each pair: (1) there is a unique last common ancestor, namely: not the common ancestor of any other pair of species in the test; (2) one species is hypoxia-tolerant and the other is hypoxia-sensitive. It should be noted that two species from the “additional group” (shrew and elephant) were used in the balanced experimental design to represent hypoxia-sensitive species. This is because these species are the most closely related species (for which p53 sequence was published) to the stare-nosed-mole and manatee, respectively, and despite some stress-related inconclusive data (for example, it is not clear whether shrew faces metabolic stress due to its small body mass).
Protein and mRNA sequences were retrieved from NCBI protein and nucleotide databases, respectively [http://www.ncbi.nlm.nih.gov/]. Habitat and behavioral information were obtained from the Animal Diversity Web (ADW) online database . Taxonomy information was collected from the NCBI Taxonomy browser . Maximal lifespan and adult weight information were obtained from AnAge database . Accession numbers, taxonomy information, lifespan, weight, and additional information related to stress adaptations are available in Additional file 1: Table S1.
Analysis of the p53 protein sequence and its functional domains
Forty-seven p53 sequences of species included in the hypoxic-, metabolic- and non-stress groups were aligned using MAFFT (multiple alignment using fast Fourier transform) L-INS-i (v7), which applies an iterative refinement method for local pairwise alignment, based on the BLOSUM62 model . To obtain p53 fragments for each functional domain in p53, the alignment result was trimmed using AliView alignment editor  to the following domains: TAD1, TAD2, PRD, DBD, NLS, and TD/RD, according to relative positions in human p53 (Fig. 1b). Gaps were removed from the sub-alignments, and the resulting fragments were used for further analysis as described below. The fragments of each p53 functional domain were aligned using MAFFT (as described above). These alignments were used to create a distance matrix according to a Poisson correction model , using MEGA6 .
Domain specific intra-group similarity
In order to test whether pairwise distances within a group involving hypoxia-tolerant species and a matched group of hypoxia-sensitive species are similar, the following test was performed. We characterized domain specific intra-group similarity for group G of k taxa by a vector S G of n = k(k-1)/2 distances, scored based on the selected domain for all pairs of taxa from this group. We then tested the following hypothesis H0: for the selected domain, components of vectors S HT and S HS have the same intra-group joint distribution. Testing this hypothesis by a standard one-tailed Mann–Whitney U test is problematic as the components of each of these two vectors are not independent. For example, information on d ij and d jk provides some information on d ik . To overcome this difficulty, we extended the standard U test by taking into account this type of dependency in row data as described below. Analogously to a standard U test, we used test statistic U equal to the number of pairs (out of all n HT *n HS ) of components of vectors S HT and S HS such that S HT > S HS . To conduct the test on real data, we estimate the distribution of U values under hypothesis H0 by the following Monte-Carlo  simulations: In each of N = 100,000 runs we simulated k HT + k HS random points in Euclidean space. Coordinates of these points were simulated as realization of independent random values with uniform distribution from 0 to 1. For groups of the first k HT and of the next k HS points, we calculated all n HT and n HS intra-group (dependent) distances as well as U-value. We further estimate p-value for the observed real vectors S HT and S HS as the proportion of simulated runs (out of N = 100,000) with U-E U ≥ U obs-E U, where E U = (n HT *n HS )/2 is an expectation of statistic U under hypothesis H0. If this p-value is less than some selected cutoff (e.g., 0.05) then the hypothesis H0 should be rejected.
Comparing Spalax p53 domains to the domains of other species
Spalax p53 domains were compared to the p53 domains of a range of species: (1) species included in the entire NCBI non-redundant (nr) database (2) selected 46 species, as described below. We used two types of similarity measures (i) Spalax p53 domains were used as Blastp queries against the nr database, which includes most of the non-redundant protein sequences. Blastp hits were then sorted by E-value; (ii) The distances of p53 domains belonging to selected 46 species, from Spalax were calculated using Poisson correction model. Fisher’s exact was used to test whether the domain-distance categories (namely, ‘close’ and ‘distant’ to Spalax) are significantly associated with the stress-tolerance categories. The categories ‘distant’ and ‘close’ were measured by: (a) expected distance by phylogeny; (b) observed domain sequence distance. An ‘expected’ 2 × 2 contingency table was built with the counts of ‘hypoxia-tolerant’/‘hypoxia-sensitive’ species that are ‘close’/‘distant’ from Spalax, where ‘close’ species were classified as those belonging to the Rodentia or Lagomorpha clades (with less than 87 million years divergence time from Spalax). In this test, ‘hypoxia-tolerant’ refers to 15 species included in the hypoxic-stress group and ‘hypoxia-sensitive’ refers to 32 species included in the metabolic- and non-stress groups (see Additional file 1: Table S1). The ‘observed’ 2 × 2 contingency table included the same total counts of ‘distant’/‘close’ species as in the ‘expected’ table, but the distance measure was domain distance rather than phylogenetic distance. Fisher’s exact test was conducted using the R function ‘fisher.test’. The data set tested included a higher proportion of species categorized as phylogenetically “distant” than phylogenetically “close” among the HT group (15 species), and a lower proportion for the HS group (32 species). Thus, Fisher’s exact test odds ratio (OR) parameter was adjusted and defined as (Hclose/Hdist)/(Sclose/Sdist); where: Hclose and Hdist are the counts of HT species that are either closely-related or distantly-related to Spalax, respectively; Sclose and Sdist are the counts of HS species that are either closely related or distantly related to Spalax, respectively. The parameter ‘or’ (odds ratio) of the R function ‘fisher.test’ was calculated from the ‘expected’ table, to reflect the expected proportions of HT species in the ‘close’ vs. ‘distant’ group (see Additional file 3: Table S2, TAD2 example). In addition, the test was repeated for a subset of 43 species and of 41 representative species (for groups of species separated by <10 million years, only a single species was selected, and the rest were excluded). Furthermore, the test was repeated 100 times with random sub-samples of 38 and 33 species, out of the 43 species (jackknife resampling method) and the mean p-value was calculated. In the above-described test we used a data set that included all known hypoxia-tolerant species for which p53 sequence was published. This comprises species within clades that consist of only hypoxia-tolerant species (Additional file 2: Figure S1), and thus it may introduce potential biases due to dependent samples. Therefore, the results were tested also in a balanced data set (as explained in the previous Methods sections).
For classification analysis, we used the data-mining method “Boosting Trees for Classification” (http://www.statsoft.com/Textbook/Boosting-Trees-Regression-Classification) . This method predicts the membership of objects in different groups based on the values of one or more predictor variables. The p53 protein sequences of 66 species (see Additional file 1: Table S1) were aligned using CLUSTALW2 Omega (1.2.1) , where each alignment position was considered as a potential predictor variable in the analysis. The actual predictors were set to the positions that were found to be more frequent in species adapted to hypoxic/metabolic stress (see Results, Fig. 2a). The species’ classification into three categories (hypoxic-, metabolic-, and non-stress) was set to be the dependent variable, and the algorithm classified the species based on the defined predictors. The “Training” data was set to be the 47 species that we subjectively assigned to the “hypoxic-stress”, “metabolic-stress”, and “non-stress” groups. This classification was determined according to information we collected from the literature related to the ability of the species to endure different types of hypoxic/metabolic stresses (see Additional file 1: Table S1). Thus, these species were used by the boosted-tree-algorithm to build the model, and their subjective classification was set to be the initial observed values. The remaining 19 species (“Additional” group) comprised the “Test” data that was used to evaluate the fitting of the model over successive iterations and did not participate in the building of the model. The initial observed values of these species were tentatively set to the intermediate stress type, namely: “metabolic stress”.
Structural analysis of Spalax p53
PyMOL Molecular Graphics System, Version 184.108.40.206 Schrödinger, LLC was used for molecular visualization of proteins’ interaction with p53. For structures involving p53 TAD2 (Additional file 9: Figure S5), all residues that differ between human and Spalax p53 sequences were marked and displayed in balls. For visualizing ‘Selecton’ positive-selection results (Fig. 4), the residues in TAD2 and RD for which dN/dS >1 were marked by the Selecton software in the relevant PDB file and then displayed by PyMOL, accordingly. For modeling the structure of Spalax p53 in complex with RPA70 according to the human complex, only the residues participating in the direct interface between p53 and RPA70 were replaced by Spalax residues and displayed. For modeling the structure of Spalax p53 in complex with p300, only residues that differ between Spalax and humans and participate in the stress-related patterns were replaced and displayed. PDB files were retrieved from the PDB database . MutaBind Server  was used to estimate the putative changes in the binding affinity of Spalax p53 to other DDR proteins, compared to other species.
Site-specific positive-selection analysis
A phylogenetic tree of 65 species was built according to the phylogenetic tree topologies in TimeTree of Life , a resource for estimations of divergence-time between species. Then, the tree-branch lengths were optimized using Rate4Site  according to the p53 protein sequences of the 65 species, aligned by MAFFT (Additional file 11: Table S6). Since site-specific positive-selection analysis tools are extremely sensitive to alignment errors, when aligning p53 codon sequences, we manually inspected, using AliView alignment editor , the results of two codon-alignment methods: (1) PAGAN v.0.47  and (2) the method implemented by the Selecton server , which is based on MUSCLE multiple alignment . Due to an ambiguous result for amino acid residues in the disordered TAD2 when using PAGAN alignment, the aligned codon sequences produced by Selecton were chosen for use as input data in the analytical procedures described below, in addition to the phylogenetic tree of 65 species for which the branches were optimized.
Positive-selection analysis was conducted using: (1) The HyPhy module MEME . This method reports site-specific dN/dS (i.e., non-synonymous vs. synonymous substitution rate), allowing the estimated dN/dS rate to change between both sites and branches, as expected in the case of p53, where it is likely that at different sites, only some branches were positively selected. For this analysis, following the HyPhy manual: SLAC analysis was run first, followed by MEME analysis, both analyses with HyPhy QuickSelectionDetection.bf module; (2) Maximum likelihood estimation of site-specific dN/dS using CODEML program in PAML version 4.8 [96, 97]. This program allows comparing the models M1a vs. M2a, and M7 vs. M8. The null-models M1a and M7 estimate the probabilities of 2 and 10 classes of dN/dS for each site, all with dN/dS between 0 and 1. The alternative models M2a and M8 include an additional site class with dN/dS >1 estimate. Using a likelihood ratio test (LRT), p-values were calculated from χ2 distribution. If a specific alternative model was found to be significant (LRT p-value < 0.05), then positive selection was tested for each site (p > 0.95 for the dN/dS >1 site class). Site class probability was tested using Bayes Empirical Bayes (BEB) . (3) The Selecton online server  [http://selecton.tau.ac.il/]. This server implements an empirical Bayesian method for dN/dS calculation and allows to compare the M8 model (that assumes selection) vs. the M8a null model (that assumes no selection) by LRT as described above, and to identify site-specific selection, accordingly. It is important to clarify that site-specific positive selection models do not provide information regarding the lineage in which the positive selection might occur. Methods that provide information regarding the site and the branch (e.g., branch-site model, PAML) were not used since these methods are not reliable for regions containing alignment gaps, which were tested here .
Nevo E, Ivanitskaya E, Beiles A. Adaptive radiation of blind subterranean mole rats: naming and revisiting the four sibling species of the Spalax ehrenbergi superspecies in Israel: Spalax galili (2n = 52), S. golani (2n = 54), S. carmeli (2n = 58) and S. judaei (2n = 60). Leiden: Backhuys Publishers; 2001.
Avivi A, Resnick MB, Nevo E, Joel A, Levy AP. Adaptive hypoxic tolerance in the subterranean mole rat Spalax ehrenbergi: the role of vascular endothelial growth factor. FEBS Lett. 1999;452(3):133–40.
Avivi A, Shams I, Joel A, Lache O, Levy AP, Nevo E. Increased blood vessel density provides the mole rat physiological tolerance to its hypoxic subterranean habitat. FASEB J. 2005;19(10):1314–6.
Shams I, Avivi A, Nevo E. Oxygen and carbon dioxide fluctuations in burrows of subterranean blind mole rats indicate tolerance to hypoxic-hypercapnic stresses. Comp Biochem Physiol A Mol Integr Physiol. 2005;142(3):376–82.
Lanoix D, Lacasse AA, Reiter RJ, Vaillancourt C. Melatonin: the watchdog of villous trophoblast homeostasis against hypoxia/reoxygenation-induced oxidative stress and apoptosis. Mol Cell Endocrinol. 2013;381(1-2):35–45.
Aguilera A, Garcia-Muse T. Causes of genome instability. Annu Rev Genet. 2013;47:1–32.
Bristow RG, Hill RP. Hypoxia and metabolism. Hypoxia, DNA repair and genetic instability. Nat Rev Cancer. 2008;8(3):180–92.
Dizdaroglu M. Oxidatively induced DNA damage: mechanisms, repair and disease. Cancer Lett. 2012;327(1-2):26–47.
Klein TJ, Glazer PM. The tumor microenvironment and DNA repair. Semin Radiat Oncol. 2010;20(4):282–7.
Pires IM, Bencokova Z, Milani M, Folkes LK, Li JL, Stratford MR, Harris AL, Hammond EM. Effects of acute versus chronic hypoxia on DNA damage responses and genomic instability. Cancer Res. 2010;70(3):925–35.
Maslov AY, Vijg J. Genome instability, cancer and aging. Biochim Biophys Acta. 2009;1790(10):963–9.
Tacutu R, Craig T, Budovsky A, Wuttke D, Lehmann G, Taranukha D, Costa J, Fraifeld VE, de Magalhaes JP. Human ageing genomic resources: integrated databases and tools for the biology and genetics of ageing. Nucleic Acids Res. 2013;41(Database issue):D1027–1033.
Manov I, Hirsh M, Iancu TC, Malik A, Sotnichenko N, Band M, Avivi A, Shams I. Pronounced cancer resistance in a subterranean rodent, the blind mole-rat, Spalax: in vivo and in vitro evidence. BMC Biol. 2013;11:91.
Buffenstein R. The naked mole-rat: a new long-living model for human aging research. J Gerontol Ser A Biol Med Sci. 2005;60(11):1369–77.
Delaney MA, Ward JM, Walsh TF, Chinnadurai SK, Kerns K, Kinsel MJ, Treuting PM. Initial case reports of cancer in naked mole-rats (heterocephalus glaber). Vet Pathol. 2016;53(3):691–6.
Hedges SB, Dudley J, Kumar S. TimeTree: a public knowledge-base of divergence times among organisms. Bioinformatics. 2006;22(23):2971–2.
Larson J, Park TJ. Extreme hypoxia tolerance of naked mole-rat brain. Neuroreport. 2009;20(18):1634–7.
de Magalhães JP. How ageing processes influence cancer. Nat Rev Cancer. 2013;13(5):357–65.
Wilhelm Filho D, Sell F, Ribeiro L, Ghislandi M, Carrasquedo F, Fraga CG, Wallauer JP, Simoes-Lopes PC, Uhart MM. Comparison between the antioxidant status of terrestrial and diving mammals. Comp Biochem Physiol A Mol Integr Physiol. 2002;133(3):885–92.
Carey HV, Frank CL, Seifert JP. Hibernation induces oxidative stress and activation of NK-kappaB in ground squirrel intestine. J Comp Physiol B. 2000;170(7):551–9.
Carey HV, Andrews MT, Martin SL. Mammalian hibernation: cellular and molecular responses to depressed metabolism and low temperature. Physiol Rev. 2003;83(4):1153–81.
Orr AL, Lohse LA, Drew KL, Hermes-Lima M. Physiological oxidative stress after arousal from hibernation in Arctic ground squirrel. Comp Biochem Physiol A Mol Integr Physiol. 2009;153(2):213–21.
Toien O, Drew KL, Chao ML, Rice ME. Ascorbate dynamics and oxygen consumption during arousal from hibernation in Arctic ground squirrels. Am J Physiol Regul Integr Comp Physiol. 2001;281(2):R572–583.
Merkt JR, Taylor CR. “Metabolic switch” for desert survival. Proc Natl Acad Sci U S A. 1994;91(25):12313–6.
el Hachimi Z, Tijane M, Boissonnet G, Benjouad A, Desmadril M, Yon JM. Regulation of the skeletal muscle metabolism during hibernation of Jaculus orientalis. Comp Biochem Physiol B. 1990;96(3):457–9.
Schulke S, Dreidax D, Malik A, Burmester T, Nevo E, Band M, Avivi A, Hankeln T. Living with stress: regulation of antioxidant defense genes in the subterranean, hypoxia-tolerant mole rat, Spalax. Gene. 2012;500(2):199–206.
Vazquez-Medina JP, Zenteno-Savin T, Elsner R, Ortiz RM. Coping with physiological oxidative stress: a review of antioxidant strategies in seals. J Comp Physiol B. 2012;182(6):741–50.
Ali MA, Kazzam E, Amir N, Nyberg F, Adem A. Effects of dehydration and blockade of angiotensin II AT1 receptor on stress hormones and anti-oxidants in the one-humped camel. BMC Vet Res. 2013;9:232.
Seim I, Fang X, Xiong Z, Lobanov AV, Huang Z, Ma S, Feng Y, Turanov AA, Zhu Y, Lenz TL, et al. Genome analysis reveals insights into physiology and longevity of the Brandt’s bat Myotis brandtii. Nat Commun. 2013;4:2212.
Maina JN. What it takes to fly: the structural and functional respiratory refinements in birds and bats. J Exp Biol. 2000;203(Pt 20):3045–64.
Salmon AB, Leonard S, Masamsetti V, Pierce A, Podlutsky AJ, Podlutskaya N, Richardson A, Austad SN, Chaudhuri AR. The long lifespan of two bat species is correlated with resistance to protein oxidation and enhanced protein homeostasis. FASEB J. 2009;23(7):2317–26.
Wilkinson GS, South JM. Life history, ecology and longevity in bats. Aging Cell. 2002;1(2):124–31.
Munshi-South J, Wilkinson GS. Bats and birds: Exceptional longevity despite high metabolic rates. Ageing Res Rev. 2010;9(1):12–9.
Zhang G, Cowled C, Shi Z, Huang Z, Bishop-Lilly KA, Fang X, Wynne JW, Xiong Z, Baker ML, Zhao W, et al. Comparative analysis of bat genomes provides insight into the evolution of flight and immunity. Science. 2013;339(6118):456–60.
Reinhardt HC, Schumacher B. The p53 network: cellular and systemic DNA damage responses in aging and cancer. Trends Genet. 2012;28(3):128–36.
Polo SE, Jackson SP. Dynamics of DNA damage response proteins at DNA breaks: a focus on protein modifications. Genes Dev. 2011;25(5):409–33.
Lim YP, Lim TT, Chan YL, Song AC, Yeo BH, Vojtesek B, Coomber D, Rajagopal G, Lane D. The p53 knowledgebase: an integrated information resource for p53 research. Oncogene. 2007;26(11):1517–21.
Rodier F, Campisi J, Bhaumik D. Two faces of p53: aging and tumor suppression. Nucleic Acids Res. 2007;35(22):7475–84.
Avivi A, Ashur-Fabian O, Joel A, Trakhtenbrot L, Adamsky K, Goldstein I, Amariglio N, Rechavi G, Nevo E. P53 in blind subterranean mole rats--loss-of-function versus gain-of-function activities on newly cloned Spalax target genes. Oncogene. 2007;26(17):2507–12.
Band M, Ashur-Fabian O, Avivi A. The expression of p53-target genes in the hypoxia-tolerant subterranean mole-rat is hypoxia-dependent and similar to expression patterns in solid tumors. Cell Cycle. 2010;9(16):3347–52.
Uversky VN, Dave V, Iakoucheva LM, Malaney P, Metallo SJ, Pathak RR, Joerger AC. Pathological unfoldomics of uncontrolled chaos: intrinsically disordered proteins and human diseases. Chem Rev. 2014;114(13):6844–79.
Xue B, Brown CJ, Dunker AK, Uversky VN. Intrinsically disordered regions of p53 family are highly diversified in evolution. Biochim Biophys Acta. 2013;1834(4):725–38.
Zou Y, Liu Y, Wu X, Shell SM. Functions of human replication protein A (RPA): from DNA replication to DNA damage and stress responses. J Cell Physiol. 2006;208(2):267–73.
Luo J, Nikolaev AY, Imai S, Chen D, Su F, Shiloh A, Guarente L, Gu W. Negative control of p53 by Sir2alpha promotes cell survival under stress. Cell. 2001;107(2):137–48.
Zhou XY, Shibusawa N, Naik K, Porras D, Temple K, Ou H, Kaihara K, Roe MW, Brady MJ, Wondisford FE. Insulin regulation of hepatic gluconeogenesis through phosphorylation of CREB-binding protein. Nat Med. 2004;10(6):633–7.
Mann HB, Whitney DR. On a test of whether one of two random variables is stochastically larger than the other. The annals of mathematical statistics. 1947;18(1):50–60.
Mas C, Lussier-Price M, Soni S, Morse T, Arseneault G, Di Lello P, Lafrance-Vanasse J, Bieker JJ, Omichinski JG. Structural and functional characterization of an atypical activation domain in erythroid Kruppel-like factor (EKLF). Proc Natl Acad Sci U S A. 2011;108(26):10484–9.
Follis AV, Llambi F, Ou L, Baran K, Green DR, Kriwacki RW. The DNA-binding domain mediates both nuclear and cytosolic functions of p53. Nat Struct Mol Biol. 2014;21(6):535–43.
Walker DR, Bond JP, Tarone RE, Harris CC, Makalowski W, Boguski MS, Greenblatt MS. Evolutionary conservation and somatic mutation hotspot maps of p53: correlation with p53 protein structural and functional features. Oncogene. 1999;18(1):211–8.
Abegglen LM, Caulin AF, Chan A, Lee K, Robinson R, Campbell MS, Kiso WK, Schmitt DL, Waddell PJ, Bhaskara S, et al. Potential Mechanisms for Cancer Resistance in Elephants and Comparative Cellular Response to DNA Damage in Humans. JAMA. 2015;314(17):1850–60.
Jenkins LM, Yamaguchi H, Hayashi R, Cherry S, Tropea JE, Miller M, Wlodawer A, Appella E, Mazur SJ. Two distinct motifs within the p53 transactivation domain bind to the Taz2 domain of p300 and are differentially affected by phosphorylation. Biochemistry. 2009;48(6):1244–55.
Shan B, Li DW, Bruschweiler-Li L, Bruschweiler R. Competitive binding between dynamic p53 transactivation subdomains to human MDM2 protein: implications for regulating the p53.MDM2/MDMX interaction. J Biol Chem. 2012;287(36):30376–84.
Edwards RJ, Palopoli N. Computational prediction of short linear motifs from protein sequences. Methods Mol Biol. 2015;1268:89–141.
Van Roey K, Uyar B, Weatheritt RJ, Dinkel H, Seiler M, Budd A, Gibson TJ, Davey NE. Short linear motifs: ubiquitous and functionally diverse protein interaction modules directing cell regulation. Chem Rev. 2014;114(13):6733–78.
Dinkel H, Michael S, Weatheritt RJ, Davey NE, Van Roey K, Altenberg B, Toedt G, Uyar B, Seiler M, Budd A, et al. ELM--the database of eukaryotic linear motifs. Nucleic Acids Res. 2012;40(Database issue):D242–251.
Davey NE, Cyert MS, Moses AM. Short linear motifs - ex nihilo evolution of protein regulation. Cell Commun Signal. 2015;13:43.
Di Lello P, Jenkins LM, Jones TN, Nguyen BD, Hara T, Yamaguchi H, Dikeakos JD, Appella E, Legault P, Omichinski JG. Structure of the Tfb1/p53 complex: Insights into the interaction between the p62/Tfb1 subunit of TFIIH and the activation domain of p53. Mol Cell. 2006;22(6):731–40.
Bochkareva E, Kaustov L, Ayed A, Yi GS, Lu Y, Pineda-Lucena A, Liao JC, Okorokov AL, Milner J, Arrowsmith CH, et al. Single-stranded DNA mimicry in the p53 transactivation domain interaction with replication protein A. Proc Natl Acad Sci U S A. 2005;102(43):15412–7.
Rowell JP, Simpson KL, Stott K, Watson M, Thomas JO. HMGB1-facilitated p53 DNA binding occurs via HMG-Box/p53 transactivation domain interaction, regulated by the acidic tail. Structure. 2012;20(12):2014–24.
Okuda M, Nishimura Y. Extended string binding mode of the phosphorylated transactivation domain of tumor suppressor p53. J Am Chem Soc. 2014;136(40):14143–52.
Miller Jenkins LM, Feng H, Durell SR, Tagad HD, Mazur SJ, Tropea JE, Bai Y, Appella E. Characterization of the p300 Taz2-p53 TAD2 complex and comparison with the p300 Taz2-p53 TAD1 complex. Biochemistry. 2015;54(11):2001–10.
Rappaport N, Twik M, Nativ N, Stelzer G, Bahir I, Stein TI, Safran M, Lancet D. MalaCards: A Comprehensive Automatically-Mined Database of Human Diseases. Curr Protoc Bioinformatics. 2014;47:1 24. 21-21 24 19.
Abramova NA, Russell J, Botchan M, Li R. Interaction between replication protein A and p53 is disrupted after UV damage in a DNA repair-dependent manner. Proc Natl Acad Sci U S A. 1997;94(14):7186–91.
Leiter LM, Chen J, Marathe T, Tanaka M, Dutta A. Loss of transactivation and transrepression function, and not RPA binding, alters growth suppression by p53. Oncogene. 1996;12(12):2661–8.
Romanova LY, Willers H, Blagosklonny MV, Powell SN. The interaction of p53 with replication protein A mediates suppression of homologous recombination. Oncogene. 2004;23(56):9025–33.
Brandsma I, Gent DC. Pathway choice in DNA double strand break repair: observations of a balancing act. Genome integrity. 2012;3(1):9.
Rajagopalan S, Andreeva A, Rutherford TJ, Fersht AR. Mapping the physical and functional interactions between the tumor suppressors p53 and BRCA2. Proc Natl Acad Sci U S A. 2010;107(19):8587–92.
Grossman SR. p300/CBP/p53 interaction and regulation of the p53 response. Eur J Biochem. 2001;268(10):2773–8.
Sermeus A, Michiels C. Reciprocal influence of the p53 and the hypoxic pathways. Cell Death Dis. 2011;2:e164.
Shams I, Avivi A, Nevo E. Hypoxic stress tolerance of the blind subterranean mole rat: expression of erythropoietin and hypoxia-inducible factor 1 alpha. Proc Natl Acad Sci U S A. 2004;101(26):9698–703.
Puca R, Nardinocchi L, Sacchi A, Rechavi G, Givol D, D’Orazi G. HIPK2 modulates p53 activity towards pro-apoptotic transcription. Mol Cancer. 2009;8:85.
Teufel DP, Freund SM, Bycroft M, Fersht AR. Four domains of p300 each bind tightly to a sequence spanning both transactivation subdomains of p53. Proc Natl Acad Sci U S A. 2007;104(17):7009–14.
Sakaguchi K, Herrera JE, Saito S, Miki T, Bustin M, Vassilev A, Anderson CW, Appella E. DNA damage activates p53 through a phosphorylation-acetylation cascade. Genes Dev. 1998;12(18):2831–41.
Liu L, Scolnick DM, Trievel RC, Zhang HB, Marmorstein R, Halazonetis TD, Berger SL. p53 sites acetylated in vitro by PCAF and p300 are acetylated in vivo in response to DNA damage. Mol Cell Biol. 1999;19(2):1202–9.
North BJ, Verdin E. Sirtuins: Sir2-related NAD-dependent protein deacetylases. Genome Biol. 2004;5(5):224.
Guarente L. Introduction: sirtuins in aging and diseases. Methods Mol Biol. 2013;1077:3–10.
Howitz KT, Bitterman KJ, Cohen HY, Lamming DW, Lavu S, Wood JG, Zipkin RE, Chung P, Kisielewski A, Zhang LL, et al. Small molecule activators of sirtuins extend Saccharomyces cerevisiae lifespan. Nature. 2003;425(6954):191–6.
Ramirez JM, Folkow LP, Blix AS. Hypoxia tolerance in mammals and birds: from the wilderness to the clinic. Annu Rev Physiol. 2007;69:113–43.
Williams MF. Morphological evidence of marine adaptations in human kidneys. Med Hypotheses. 2006;66(2):247–57.
Myers P, Espinosa R, Parr CS, Jones T, Hammond GS, Dewey TA. The Animal Diversity Web (online). 2015. Accessed at http://animaldiversity.org.
Sayers EW, Barrett T, Benson DA, Bryant SH, Canese K, Chetvernin V, Church DM, DiCuccio M, Edgar R, Federhen S, et al. Database resources of the National Center for Biotechnology Information. Nucleic Acids Res. 2009;37(Database issue):D5–15.
Katoh K, Standley DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013;30(4):772–80.
Larsson A. AliView: a fast and lightweight alignment viewer and editor for large datasets. Bioinformatics. 2014;30(22):3276–8.
Zuckerkandl E, Pauling L. Evolutionary divergence and convergence in proteins. Evolving genes and proteins. 1965;97:97–166.
Tamura K, Stecher G, Peterson D, Filipski A, Kumar S. MEGA6: Molecular Evolutionary Genetics Analysis version 6.0. Mol Biol Evol. 2013;30(12):2725–9.
Fishman G. Monte Carlo: concepts, algorithms, and applications. New York: Springer; 1996.
Friedman J, Hastie T, Tibshirani R. The elements of statistical learning, vol. 1. Berlin: Springer series in statistics Springer; 2001.
Sievers F, Wilm A, Dineen D, Gibson TJ, Karplus K, Li W, Lopez R, McWilliam H, Remmert M, Soding J, et al. Fast, scalable generation of high-quality protein multiple sequence alignments using Clustal Omega. Mol Syst Biol. 2011;7:539.
Bernstein FC, Koetzle TF, Williams GJ, Meyer Jr EF, Brice MD, Rodgers JR, Kennard O, Shimanouchi T, Tasumi M. The Protein Data Bank: a computer-based archival file for macromolecular structures. J Mol Biol. 1977;112(3):535–42.
Li M, Simonetti FL, Goncearenco A, Panchenko AR. MutaBind estimates and interprets the effects of sequence variants on protein-protein interactions. Nucleic Acids Res. 2016;44(W1):W494–501.
Pupko T, Bell RE, Mayrose I, Glaser F, Ben-Tal N. Rate4Site: an algorithmic tool for the identification of functional regions in proteins by surface mapping of evolutionary determinants within their homologues. Bioinformatics. 2002;18 Suppl 1:S71–77.
Loytynoja A, Vilella AJ, Goldman N. Accurate extension of multiple sequence alignments using a phylogeny-aware graph algorithm. Bioinformatics. 2012;28(13):1684–91.
Stern A, Doron-Faigenboim A, Erez E, Martz E, Bacharach E, Pupko T. Selecton 2007: advanced models for detecting positive and purifying selection using a Bayesian inference approach. Nucleic Acids Res. 2007;35(Web Server issue):W506–511.
Edgar RC. MUSCLE: a multiple sequence alignment method with reduced time and space complexity. BMC bioinformatics. 2004;5:113.
Murrell B, Wertheim JO, Moola S, Weighill T, Scheffler K, Kosakovsky Pond SL. Detecting individual sites subject to episodic diversifying selection. PLoS Genet. 2012;8(7):e1002764.
Yang Z. PAML 4: phylogenetic analysis by maximum likelihood. Mol Biol Evol. 2007;24(8):1586–91.
Yang Z, Nielsen R. Mutation-selection models of codon substitution and their use to estimate selective strengths on codon usage. Mol Biol Evol. 2008;25(3):568–79.
Yang Z, Wong WS, Nielsen R. Bayes empirical bayes inference of amino acid sites under positive selection. Mol Biol Evol. 2005;22(4):1107–18.
Fletcher W, Yang Z. The effect of insertions, deletions, and alignment errors on the branch-site test of positive selection. Mol Biol Evol. 2010;27(10):2257–67.
We wish to thank Adi Stern and Kazutaka Katoh for their kind help in using Selecton and MAFFT servers, respectively, Svetlana Frenkel for her help in allocating p53 mRNA sequences, and Robert Altwasser for his valuable advices. We wish to give special acknowledgment to Professor David Givol (from the Weizmann Institute of Science, Israel), of blessed memory, for thought-provoking discussions and scientific support in the initial stage of the study.
This work was supported by John Templeton Foundation [grant #5305] and by the generous donation of the Kadas Family (UK).
Availability of data and material
All data generated or analyzed during this study are included in this published article and its supplementary information files.
Wrote the paper: VD, AK, IS, AM. Analyzed the data: VD, AM, AK, ZF. Performed structural analyses: YO. Designed the experiments: AK, VD, AM, IS, AA, IM. Collected the data: VD. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
Additional information on the species employed in the study. (XLS 56 kb)
Phylogenetic tree of 47 species. (PDF 15 kb)
P53 sequences distances. (XLS 90 kb)
Stress related patterns and classification analysis. (XLS 83 kb)
Stress related patterns in p53-TAD2 and KLF1-TAD2. (PDF 23 kb)
Phylogenetic tree of 66 species. (PDF 16 kb)
Classification analysis. (PDF 4479 kb)
Critical DDR proteins that interact with TAD2 and/or TD/RD. (XLS 42 kb)
Differences between human and Spalax in p53-TAD2 interactions. (PDF 1461 kb)
Positive selection results. (XLS 443 kb)
Positive selection input. (XLS 269 kb)
About this article
Cite this article
Domankevich, V., Opatowsky, Y., Malik, A. et al. Adaptive patterns in the p53 protein sequence of the hypoxia- and cancer-tolerant blind mole rat Spalax . BMC Evol Biol 16, 177 (2016). https://doi.org/10.1186/s12862-016-0743-8
- Oxidative stress
- Convergent evolution
- Cancer resistance