32125 Kinetic profiling of protein-drug binding using non-equilibrium molecular dynamics simulations
Begeleider(s): Wouter Vervust

Richtingen: Master of Science in Biomedical Engineering, Master of Science in Engineering Physics

Probleemstelling:

Targeted drug design: thermodynamic and kinetic selectivity

Targeted drug design involves creating drug molecules that selectively bind to a specific target protein. Historically, the design process has relied on identifying ligands that fit into the target protein's binding pocket and have a high binding strength, as measured by the dissociation constant Kd. However, this approach has limitations because it does not consider the free energy barrier that must be overcome during the binding reaction, and it neglects the dynamics along the reaction coordinate. As a result, drug molecules that have a high affinity for the target protein can still bind to off-target proteins or have poor overall selectivity. To overcome these limitations, the concept of kinetic selectivity, which considers the rate constants kon and koff, has become increasingly important in recent years. By focusing on both thermodynamic and kinetic selectivity, drug designers can more effectively develop compounds that bind selectively and with high residence time to their intended targets, minimizing off-target toxicity and improving overall efficacy. A schematic picture of protein-ligand binding is given in Figure 1.

Figure 1: Left: schematic picture of the protein-ligand (P-L) binding free energy profile. The equilibrium constant Kd is fully described by the free energy difference ΔG of the unbound and bound states, describing the strength of the P-L bound complex. For this reason, Kd has served as the dominant parameter in targeted drug design. However, the dynamics of the binding process are often not negligible, and can highly impact the effectiveness of the candidate ligand in vivo (see text). Right: rate equation accompanying the bimolecular binding reaction, where [P] and [L] denote the protein and ligand concentrations, respectively.

Molecular dynamics for protein-ligand binding

Computational techniques have greatly advanced high-throughput applications, such as protein-ligand docking in targeted drug design, and molecular dynamics (MD) simulations have enabled the fully atomistic depiction of proteins. However, the binding/residence times are typically multiple milliseconds, seconds, minutes, or even hours, while conventional MD simulations rarely exceed one millisecond. Thus, binding interactions, which are rare and/or slow events in simulations, pose a challenge in studying time-dependent parameters such as the (un)binding kinetics. This thesis aims to develop advanced methodological simulations that artificially accelerate the binding process to retrieve the true dynamics of protein-ligand binding using principles from statistical physics.

Recently, two methods have been developed to extract both the thermodynamics and kinetics of rare events: the Transition-based Reweighting Analysis Method (TRAM) [1] and the Dynamic Histogram Analysis Method Extended to Detailed balance (DHAMed) [2]. Both methods employ a Markov State Model, where the simulation trajectory data is discretized in predefined phase space states (TRAM) or along a predefined reaction coordinate (DHAMed). Compared to earlier methods, both TRAM and DHAMed improve Markovianity, which assumes that each transition between states/coordinates depends only on the state/coordinate it currently occupies, and not on its history.

However, due to the need for low-dimensional reaction coordinates, simulations struggle to capture binding pathways that contain energy barriers in degrees of freedom orthogonal to the reaction coordinate. This creates a problem with the Markov assumption, as transitions become dependent on their binding pathway rather than solely on their position along the reaction coordinate. To address this, both TRAM and DHAMed employ the Markov model by incorporating "local" equilibrium sampling instead of "global" equilibrium sampling used in older methods like WHAM and MBAR. Overall, by employing advanced methodological simulations and statistical physics principles, this thesis aims to derive the kinetic rate constants of protein-ligand binding and provide a more accurate understanding of the binding process, which will be valuable for drug design and related fields.


Doelstelling:

Objectives

This thesis will use the Abelson tyrosine kinase (ABL) as a test case to study the role of orthogonal DoFs (Figure 2). This protein is involved in cell differentiation, cell division, cell adhesion and stress response [3]⁠. Some patients develop the BCR-ABL oncoprotein which is always active, resulting in chronic myeloid leukaemia [4]⁠. Small drug molecules (e.g. imatinib) have been developed that force BCR-ABL in the inactive state, by occupying the ATP-binding site of ABL.

The student will get familiar with the ABL-imatinib system, as well as the theoretical concepts of biased simulations and Markov state models. The student will then apply the DHAMed/TRAM methodology using a straight-forward reaction coordinate, after which the latter can be optimized. Depending on the creative interest of the student, a new methodology tailored for the treatment of the orthogonal DoFs in binding kinetics can be developed. In short, the student will run and analyze advanced molecular dynamics simulations on the ABL-ligand complex to test whether the binding kinetics can be retrieved.

Figure 2: left: imatinib (magenta) bound to ABL. Right: different conformations of ABL, where important flexible loops are colored red/orange (P-loop, active/inactive) and blue/seablue (A-loop, inactive/active)

References

[1] Wu, H., Paul, F., Wehmeyer, C., & Noé, F. (2016). Multiensemble Markov models of molecular thermodynamics and kinetics. Proceedings of the National Academy of Sciences, 113(23), E3221-E3230.

[2] L. S. Stelzl, A. Kells, E. Rosta, and G. Hummer, “Dynamic Histogram Analysis To Determine Free Energies and Rates from Biased Simulations,” J. Chem. Theory Comput., vol. 13, no. 12, pp. 6328–6342, 2017, doi: 10.1021/acs.jctc.7b00373.

[3] J. Y. J. Wang, “The Capable ABL: What Is Its Biological Function?,” Mol. Cell. Biol., vol. 34, no. 7, pp. 1188–1197, 2014, doi: 10.1128/mcb.01454-13.

[4] A. Hai, N. A. Kizilbash, S. H. H Zaidi, J. Alruwaili, and K. Shahzad, “Differences in structural elements of Bcr-Abl oncoprotein isoforms in Chronic Myelogenous Leukemia,” Bioinformation, vol. 10, no. 3, pp. 108–114, 2014, doi: 10.6026/97320630010108.