An Integrated View of p53 Dynamics, Function, and Reactivation

The tumor suppressor p53 plays a vital role in responding to cell stressors such as DNA damage, hypoxia, and tumor formation by inducing cell-cycle arrest, senescence, or apoptosis. Expression level alterations and mutational frequency implicates p53 in most human cancers. In this review, we show how both computational and experimental methods have been used to provide an integrated view of p53 dynamics, function, and reactivation potential. We argue that p53 serves as an exceptional case study for developing methods in modeling intrinsically disordered proteins. We describe how these methods can be leveraged to improve p53 reactivation molecule design and other novel therapeutic modalities, such as PROteolysis TARgeting Chimeras (PROTACs).


The Opportunities and Challenges of the Tumor Suppressor p53
The transcription factor p53 acts as arguably the most important tumor suppressor protein in humans, with more than 50% of all human cancers implicating mutated p53 in their progression. This transcription factor regulates key genes that are involved in cell function and carcinogenesis, such as apoptosis, senescence, and DNA repair. Because of mutant p53's frequent involvement in cancer progression, it is thought of as a promising and exciting drug target in the field of cancer therapy [1,2]. However, many p53 mutants are known to have significantly altered dynamics and function in comparison to the wildtype protein, challenging traditional structure-based drug discovery methods in the quest for mutant p53 reactivation. Furthermore, p53 consists of multiple domains, adding an additional level of complexity to the characterization of this protein: the intrinsically disordered N-terminal transactivation domain (NTD, residues 1-93) that is responsible for activating target genes, the sequence-specific DNA binding domain (DBD, residues 94-294), the tetramerization domain (TET, residues 323-360) which assists in formation of the functional p53 tetramer, and the intrinsically disordered C-terminal regulatory domain (CTD, residues 363-393) ( Figure 1).
In this review, we describe how modern computational methodologies, tightly integrated with experiment, are emerging as capable of providing avenues for reactivation of p53 cancer mutants. We discuss recent computational and experimental advances in characterizing each of these domains separately as well as in the context of the larger oligomeric, quaternary complex. Lastly, we discuss the therapeutic implications in terms of mutant p53 reactivation and protein-protein interactions, and provide a perspective on how the approaches developed to study this important system can be useful to guide work on similar challenging systems, including intrinsically disordered proteins (IDPs) more generally and in the context of emerging therapeutic modalities, such as PROteolysis TARgeting Chimeras (PROTACs). This mini-review covers the most recent developments in p53 reactivation research; for more comprehensive reviews, we refer readers to Refs. [1][2][3][4][5].

Insight into individual p53 domains
Despite lacking a fixed three-dimensional structure, intrinsically disordered proteins (IDPs) can adopt one or more secondary structures upon target binding (protein, small molecule, RNA, DNA or ions), post-translational modifications or chemical environment changes. Given the inherent instability of p53, its multidomain architecture with long intrinsically disordered stretches and the challenges it poses to experimental studies, molecular dynamics (MD) simulations have provided a useful window into the function and dynamics of the wildtype p53 protein and its oncogenic mutants. Computational costs associated with allatom simulations, however, limit the size of the systems under investigation, and thus, the great majority of studies have focused on individual protein domains.

DBD
The p53 DBD is intrinsically unstable, often aggregating with small perturbations in temperature, ionic environment, or mutations. This domain is a hotspot region for cancer mutations ( Figure 1) and has been the focus of several recent computational studies highlighting the diverse effects of the single-point mutations. MD simulations of several structural mutants evidenced distinct conformational sampling in the S6/S7 loop, located at the vicinity of an aggregation-prone region, compared to wildtype [6]. Changes in the loop's interaction network due to the mutations are suggested to drive the destabilization of the DBD, leading to unfolding and protein aggregation. The presence of allostery within the DBD that explains the inactivating effect of mutations located distant from the DNA binding surface has also been proposed for the R249S mutant using MD simulations [7] and the Y220C mutant following extensive sampling and Markov state models (MSMs) [8].
The flexibility of the DBD was highlighted in MD simulations of 20 mutants, which identified exacerbation of wildtype's inherent structural vulnerabilities by the mutations and suggested a general mechanism for p53 rescue by targeting the so-called "structuraldisruption motifs" [9]. Using NMR in combination with equilibrium denaturation experiments, Bej et al [10] probed the backbone dynamics of the DBD, showing a linear correlation between the equilibrium denaturation parameters with the extent of changes in conformational entropy. Upon studying aggregation of p53 mutants using MD and coarsegrained simulations of a DBD peptide segment containing the aggregation nucleating region in conjunction with biophysical experiments, Lima et al found polymorphisms in the β-sheet peptide aggregates [11]. Distinct structural destabilization and aggregation properties were also evident in experimental and computational studies of the DNA-contact R273 mutants, indicating that conformational instability is widespread among p53 mutants and the loss of contact with DNA may not be the single molecular basis of altered p53 function [12].
Significant research aims to develop wildtype p53 reactivating small molecules through the stabilization of the native, folded state over the unfolded, aggregated state [13]. Two small molecules, PRIMA-1 and its analogue APR-246, were identified in a cell-based screen from the National Cancer Institute library [14,15]. Both molecules have been shown to target mutant p53 proteins, restoring wild-type p53 transcriptional activity, and are currently in clinical trials [16]. Another class of stabilizers, termed cavity binders, have been developed to target the elongated surface crevice created by the Y220C mutation [5,17,18]. The conformational sampling afforded by computer simulations has additionally enabled the identification of promising druggable pockets in the DBD in a number of studies [6,8,19]. This includes cryptic pockets that are absent from the experimentally resolved structures but emerge in all-atom simulations, providing new therapeutic opportunities for p53 reactivation [8,19]. (Figure 2). Among these, the L1/S3 pocket is especially important as it exists across all mutant forms of p53 presenting a broad therapeutic opportunity [19].

CTD
The p53 CTD exhibits properties of an intrinsically disordered protein, binding a number of proteins by adapting various conformations: an alpha-helix, a beta-strand, a beta-turn, or disordered structures [20][21][22][23][24] (Figure 3). Kumar et al. probed the CTD's conformational behavior upon changing protein environment, finding a distinct and predictable switch of residues 380-388 between an α-helical and β-sheet secondary structures [25]. Via circular dichroism spectroscopy and FRET, they probed the importance of temperature and hydrophobic interactions with binding partners in the secondary structure of the CTD, finding that hydrophobicity allows the CTD to preferentially retain a more α-helical structure and higher temperatures allow for acquisition of secondary structure [25].
Kannan et al. investigated the binding of the p53 CTD to 5 of its binding partners, namely S100B, cyclinA, CBP, sirtuin, and Set9 using equilibrium and non-equilibrium MD in explicit solvent [26]. The equilibrium MD demonstrated that the free p53 CTD peptide fluctuates between various conformations, including the conformations captured in crystal structures of its complexes [26]. Using non-equilibrium MD, the unbinding of the p53 CTD peptide was achieved, indicating that long-range electrostatic interactions (e.g., farther than 10 Å) results in formation of reactive conformations. This work further showed how peptide folding at or close to the binding interface steers the reactive conformations toward the binding partner [26].
Two recent separate studies on the p53 CTD peptide investigated secondary structure formation called molecular recognition features (MoRFs) with microsecond-timescale conventional MD or nanosecond-timescale replica-exchange MD [25,27]. Iida et al. focused on the interaction of the CTD peptide with S100B using virtual-system coupled multicanonical MD, a generalized ensemble MD simulation for enhanced conformational sampling [28]. They concluded that the CTD peptide adopts various conformations upon binding S100B rather than the single α-helical conformation observed in the NMR structure of the p53CTD-S100B complex. Analyzing the multimodal structural distribution of the complex, they found the conformation observed in the NMR model is the most probable orientation in the ensemble. They also found that the entropy of p53 CTD peptide in the S100B-bound state is higher than that of the free state, indicating conformational entropy may not affect p53CTD-S100B complex formation.

NTD
In the intrinsically disordered NTD, residues 10 to 40 are particularly important as they are found in a stable α-helix while in complex with MDM2, an important negative regulator of p53 [29] (Figure 3). The interface of the p53 NTD and MDM2 complex is an important drug target; as such, NTD structure and dynamics has received the attention of many researchers. Herrero-Nieto et al. combined MD at the millisecond timescale with MSMs to explore the ensemble of conformations of the p53 NTD peptide in isolation [30]. They found multiple states enriched in secondary structure elements including an α-helix as well as β-sheets in the conformational landscape of p53 NTD corresponding to about 40% of the equilibrium population. Apart from these structurally and kinetically diverse ordered states, the remaining 60% was completely heterogeneous and lacked any secondary structure. Zhao et al. coupled replica exchange molecular dynamics (REMD) simulations with MSMs to investigate the effect of dual phosphorylation at Ser46 and Thr55 on conformational ensemble of the p53 NTD peptide and found a slightly larger α-helical content in the ensemble after phosphorylation [31].
p53 and MDM2 engage in a negative autoregulatory loop, where upregulation of p53 leads to increased expression of MDM2, which in turn binds and targets p53 for ubiquitindependent degradation. Many tumors exploit this negative feedback loop by introducing MDM2 mutations to increase p53 degradation in cancer cells. Due to the co-dependence of these two proteins, there have been great strides in utilizing computational methods to understand this protein-protein interaction. For example, metadynamics simulations [32], potential of mean force (PMF) studies [33], parallel cascade selection molecular dynamics (PaCS-MD) [34], MSMs [35][36][37], and modeling employing limited data (MELD)accelerated MD [38] have been used to understand the binding/unbinding process of p53 and MDM2 at a molecular level. These methods were able to accurately predict the thermodynamic and kinetics of p53/MDM2 binding and reveal states along the p53/MDM2 binding trajectory, which can all be used as models in p53/MDM2 structure-based drug design.
Our enhanced understanding of the p53/MDM2 interaction aids in designing inhibitors to disrupt this protein-protein interaction. This is evidenced by the number of drug discovery programs targeting p53/MDM2 as discussed in a recent review by Miller, Gaiddon, and Storr [4]. The Amgen compound, AMG232 (currently in Phase I/II clinical trials), is one of the most potent p53/MDM2 inhibitors to date exhibiting sub-nanomolar binding, inducing significant p53 upregulation, and demonstrating clinically suitable PK properties. Aileron Therapeutics developed a stapled peptide, ALRN-6924, a highly potent dual inhibitor of both MDM2 and MDMX that has been shown to restore p53 activity in vitro and in vivo and is undergoing Phase I clinical studies. MDMX is also a negative regulator of p53 that binds its transactivation domain to inhibit activation of p53 transcription factors; in particular, MDMX binds MDM2 to prevent its auto-ubiquitination, which stabilizes p53 ubiquitin-dependent degradation [39,40]. Therefore, dual targeting of both MDM2 and MDMX has provided an interesting avenue in significantly enhancing anticancer activity.

Modeling the full-length p53
Structural characterization of full-length proteins that contain folded and disordered domains is a major challenge. Thus, experimental and computational approaches generally analyze various fragments of such multidomain IDPs in isolation. However, folded or disordered components of multidomain IDPs do not function as isolated entities: instead, all components of the entire protein act in synergy. In support of this, Wright et al. recently showed the NTD's importance in DNA binding accuracy using NMR and inteins, evidencing that the NTD reduces DBD binding to nonspecific p53 DNA by five-fold. This suggests competition between NTD and non-specific DNA in binding to the DBD while leaving response-element DNA binding uninhibited [41]. Furthermore, He et al. probed the NTD-DBD interactions via NMR, showing small chemical shifts suggesting weak interactions, perhaps regulating binding between the DBD and response element DNA by nucleic acid mimicry or electrostatic screening [42].
Although there are high-resolution structures available for different domains of the fulllength protein (fl-p53), there is no x-ray crystallographic or NMR structure of the entire fl-p53 [43]. Characterizing the structural ensembles and long timescale dynamics of the fl-p53 structure, or at least of large multi-domain fragments, while in complex with DNA or other protein interactors in atomic detail would be invaluable for understanding the p53-activated tumor suppression pathways, as well as how its behavior changes with p53 cancer mutations, post-translational modifications, or binding of small molecules or protein interactors. A first glimpse into the fl-p53 monomer concerted dynamics has been provided by MD simulations of Chillemi et al. [44]. These simulations revealed correlated motions between the transactivation domain 1 (TAD1) and the proline-rich region of the NTD, the tetramerization region in the CTD, and Lys120 in the DBD.
Advances in computing architectures and software performance have made possible simulations of the explicitly solvated full-length DNA-complexed p53 at sufficient timescales that can reveal domain dynamics and allosteric communication in the fl-p53 tetramer context. MD simulations of our all-atom integrative fl-p53 tetramer models bound to three different DNA sequences, namely a p21 response element, a PUMA response element and a nonspecific DNA sequence, yielded final structures consistent with electron microscopy maps and, for the first time, showed the direct interactions of the p53 CTD with DNA at the atomic level [45]. (Figure 4). Solvent mapping analyses of these nearly 1 μs-long simulations revealed multiple potential druggable pockets in p53 and a collective principal component analysis identified sequence-dependent differential quaternary binding modes of the p53 tetramer to DNA [45]. In a subsequent study including wild-type p53 and the R175H mutant in the fl-p53 tetramer context, the symmetric quaternary DNA binding mode observed for the wild-type DBDs was found absent in the R175H mutant system, similar to the case of wild-type p53 DBDs binding to nonspecific DNA in the previous study [46]. These findings prove the requirement of functional p53 and an optimal DNA sequence for a productive binding mode of DBD and DNA at the molecular level.
Methods that can identify correlated motions and allosteric communication in protein structures can be extremely helpful for the interpretation of simulations of dynamic multidomain proteins such as full-length p53. A recent development of the dynamical network analysis method by Melo et al. provides exciting opportunities for the efficient and automated identification of communication pathways in biomolecular structures, including visualization in VMD [47]. The methodology is optimized enough to efficiently process large multi-subunit protein simulations, making it exceptionally poised for p53 and other complex proteins analysis. Another method developed by Porter et al. looks for cooperative changes in solvent exposure as indication of functionally relevant conformational changes and led to the identification of not only allosteric signaling but also cryptic pockets in tested systems [48]. Furthermore, connecting dynamical network analysis with MSMs is likely to indicate how residue-level communication networks shift in different states. Multiscale methods that intersect computational approaches operating at different scales are tools capable of addressing longer-range structure, dynamics and function questions in this and other systems [49].

Conclusions & Outlook
In this review we discuss recent efforts that have provided hitherto unseen insights into the structure, dynamics, function, and druggable new pockets of p53 through integrating cutting-edge computational and experimental biophysical approaches. Techniques that drill down into one p53 domain have provide value in their own right, but structure based drug design has great potential to now also use such integrated methods on the full-length p53 structure as well as its complex with DNA. We posit that the larger-scale approaches to p53 reactivation will be particularly useful in drug design, especially in the way that they provide a more complete understanding of how dynamics affects function, including insights on longer-range allosteric effects.
Taking a higher level view, the field of computational chemistry has made significant strides in developing methods to understand protein complexes, with p53 commonly used as an appropriate example. A related field that stands to benefit from such approaches is PROTAC design. PROTACs are a promising therapeutic modality, which selectively targets proteins for degradation through exploitation of the intracellular ubiquitin-proteasome machinery [50]. These hetero-bifunctional molecules consist of a target protein ligand, an E3 ligase ligand, and an often flexible linker. This new approach in drug design shows great promise in selectively and rapidly degrading target proteins, thereby inhibiting downstream signals associated with disease pathways. Indeed, PROTAC design may be a worthwhile approach in targeting p53-associated tumors. For example, Li et al. developed MD-224, the first MDM2 PROTAC that rapidly degrades MDM2 at concentrations < 1 nM in human leukemia cells and thereby prevents p53 degradation by MDM2 [51]. Similar computational approaches, as we and others have demonstrated for p53, can be used to learn the optimal ternary complex for PROTACS, as well as understand multidomain dynamics and communication in IDPs.
Similar to p53, many multi-domain hub proteins have long IDP regions which are used for mediating the interactions with their protein partners. Impairment of such interactions are often significant and disease associated. The computational approaches used to explore fl-p53 dynamics can be adopted to investigate dynamics and long-range communication of large multi-domain proteins and their complexes. With larger and larger structural information emerging from recent advances in single particle cryoelectron microscopy and cryoelectron tomography, computational approaches that can explore such realistic complexes without losing atomic detail will become an asset.