Next Article in Journal
Parasitological and Biochemical Efficacy of the Active Ingredients of Allium sativum and Curcuma longa in Schistosoma mansoni Infected Mice
Previous Article in Journal
Intergenerational Transmission of Resistance of Callosobruchus maculatus to Essential Oil Treatment
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Molecules
Volume 26
Issue 15
10.3390/molecules26154540
molecules-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Molecular Dynamics Simulations in Designing DARPins as Phosphorylation-Specific Protein Binders of ERK2

by
Vertika Gautam
1,
Piyarat Nimmanpipug
2,3,
Sharifuddin Md Zain
1,
Noorsaadah Abd Rahman
1 and
Vannajan Sanghiran Lee
1,3,*
1
Department of Chemistry, Faculty of Science, University of Malaya, Kuala Lumpur 50603, Malaysia
2
Department of Chemistry, Faculty of Science, Chiang Mai University, Chiang Mai 50200, Thailand
3
Center of Excellence for Innovation in Analytical Science and Technology (I-ANALY-S-T), Chiang Mai University, Chiang Mai 50200, Thailand
*
Author to whom correspondence should be addressed.
Molecules 2021, 26(15), 4540; https://doi.org/10.3390/molecules26154540
Submission received: 8 July 2021 / Revised: 20 July 2021 / Accepted: 23 July 2021 / Published: 27 July 2021
(This article belongs to the Section Computational and Theoretical Chemistry)
Browse Figures
Versions Notes

Abstract

:
Extracellular signal-regulated kinases 1 and 2 (ERK1/2) play key roles in promoting cell survival and proliferation through the phosphorylation of various substrates. Remarkable antitumour activity is found in many inhibitors that act upstream of the ERK pathway. However, drug-resistant tumour cells invariably emerge after their use due to the reactivation of ERK1/2 signalling. ERK1/2 inhibitors have shown clinical efficacy as a therapeutic strategy for the treatment of tumours with mitogen-activated protein kinase (MAPK) upstream target mutations. These inhibitors may be used as a possible strategy to overcome acquired resistance to MAPK inhibitors. Here, we report a class of repeat proteins—designed ankyrin repeat protein (DARPin) macromolecules targeting ERK2 as inhibitors. The structural basis of ERK2–DARPin interactions based on molecular dynamics (MD) simulations was studied. The information was then used to predict stabilizing mutations employing a web-based algorithm, MAESTRO. To evaluate whether these design strategies were successfully deployed, we performed all-atom, explicit-solvent molecular dynamics (MD) simulations. Two mutations, Ala → Asp and Ser → Leu, were found to perform better than the original sequence (DARPin E40) based on the associated energy and key residues involved in protein-protein interaction. MD simulations and analysis of the data obtained on these mutations supported our predictions.

1. Introduction

Protein kinases play a principal regulatory role in nearly all aspects of cell biology. The human genome encodes 538 protein kinases [1]. Out of all the post-translational modifications (PTMs), protein phosphorylation is the most widespread class used in signal transduction. One of the members of the protein kinase family, mitogen-activated protein kinase (MAPK) is a type of protein kinase that is specific to the amino acids serine and threonine (i.e., a serine/threonine-specific protein kinase). MAPK is involved in the most fundamental pathway to cell biology, known as the MAPK pathway, and it plays a crucial role in integrating cell surface signals to transcriptional regulation of the proteome [2]. The MAPK pathway is also referred to as the RAS-RAF-MEK-ERK signal cascade. The main function of this cascade is to regulate physiological processes by transmitting upstream signals to its downstream effectors. They are mainly involved in cell differentiation, proliferation, survival, and death. Targeting the MAPK pathway is thought to be a hopeful strategy for cancer therapy as it is the most often mutated signalling pathway in human cancer. In the past decades, extensive efforts have been made by different research groups, leading to the clinical success of BRAF and MEK inhibitors.
In the early 1980s, the first protein kinase inhibitor was developed followed by FDA approval of many more drugs as kinase inhibitors for the treatment of cancers, for example, lung and breast cancers. Till date, around 150 drugs targeting kinases have undergone clinical phase trials, and many other specific inhibitors are in the preclinical stage of drug development [3]. However, with the speedy increase in resistance developed to clinical RAF and MEK inhibitors, interest has been encouraged towards targeting ERK directly for cancer therapy.
Natural proteins and antibodies undoubtedly changed medicine, but they usually treat only one disease target. Antibodies have some well-known limitations such as expensive production, difficult formulation, low tissue penetration, and complex architecture, and they bind to their target bivalently. On the other hand, designed ankyrin repeat proteins (DARPins) have completely overcome these limitations of conventional therapeutic approaches due to their small size, high stability, high potency, high affinity (strong binding), and flexible architecture. Moreover, the development of in vitro selection technologies such as ribosome display [4] and phage display [5] has enabled the selection of DARPins as specific binders. DARPins show commendable stability (thermodynamic as well as intracellular), which has made them promising candidates for therapeutic applications [6,7,8,9,10].
Based on naturally occurring ankyrin repeat proteins, DARPins are emerging a promising new class of binding proteins [11,12]. DARPins (designed ankyrin repeat proteins) are one of the most profusely found binding proteins in the human genome [13]. The structure of ankyrin repeat proteins consists of tightly joined repeats of 33 amino acid residues [14]. The basic structural unit of ankyrin repeats consists of two antiparallel α-helices preceding a β-turn. In a single protein, up to 29 consecutive repeats can be found [15,16]. DARPins repeat a typical structure that consists of a module flanked by N and C caps (Figure 1A), where N- and C- designate N and C terminal capping repeats, respectively, and “–” stands for the number of library modules that ranges between 2 and 4 (N2C and N3C) [17]. The hydrophobic core of “repeats” is shielded by these caps. Domains of ankyrin repeats forming a continuous hydrophobic core together with a large solvent-accessible surface that holds the repeat modules together provide stability to the structure. They are thermodynamically very stable. Through a consensus approach, it was found that repeats of DARPin consist of residues that are responsible for the maintenance of its structure (fold conservation) called fixed framework residues, while there are other residues through which DARPins interact with their target proteins known as randomized interacting residues (Figure 1B) [14]. Variations in DARPins can be brought about through randomized residues along with conserved interfaces that are present between repeat units. These interfaces are the places where single repeats can be deleted, inserted, or exchanged, maintaining the tertiary structure intact [18]. In addition to designing a large DARPin library, the consensus design approach also produced desired DARPins with enhanced properties in terms of expression levels, stability, and solubility [12,17]. All these qualities make them a suitable candidate for ERK inhibition. In the present study, we have explored DARPins (designed ankyrin proteins) as ERK2 inhibitors.
In the pharmaceutical industry, computational drug design has played a vital role in the discovery, design, and analysis of drugs. Computer technology nowadays is so rich and advanced that the accuracy of biomolecular simulations is consistently high enough to be used to truly drive preclinical drug discovery projects. Among the computational tools used for drug discovery, molecular dynamics simulations (MDS) and related methods are routinely used nowadays. Their main contribution is towards the understanding of structural flexibility together with entropic effects of complex systems. This allows an in-depth estimation of the thermodynamics and kinetics related to drug-target recognition and binding [18]. MD simulations have been extensively used in the study of protein-protein and protein-ligand interactions, and the study of the mechanism of drug action [19,20,21,22,23]. The plethora of applications of molecular dynamics simulations extend from the study of complex and dynamic processes that play a central role in biological systems to the structure determination from X-ray crystallography and NMR experiments. In biological systems, the main application of MDS focuses on studies related to the stability of proteins through conformational changes and folding. MD studies also enable the molecular recognition of cellular components such as DNA, membranes along with complexes (drug–receptor), and ion transport [24,25,26,27]. MD studies have also enabled the study of the mechanism of drug resistance [28,29,30,31]. Molecular dynamics simulations calculate the binding free energy, which is very helpful in investigating receptor-ligand interactions [32]. Over time, MDS have made immense contributions towards drug discoveries, exploring ligands such as small molecules, chemicals derived from plants, peptides, and proteins against targets such as protein kinases (PKs) [33,34,35,36,37], G-protein-coupled receptors (GPCRs) [38], and NMDA receptors [20,39,40]. Altogether, this method affords a means for drug design by providing a holistic approach to understanding the mechanism of receptor activation/deactivation, inhibiting the receptor to the mechanism of drug resistance.
During MD simulations, configurations of the progressing system are sequentially generated; these configurations come out as trajectories that contain specific details of the positions and velocities of a particle over the time of simulation. These trajectories are explored to calculate a variety of properties, such as kinetic measures and free energy and other macroscopic quantities. These properties can be further compared with experimental observables that are helpful in drug design. This method was initially perceived within theoretical physics in the late 1950s, and its application was extended to chemical physics, materials science, biomolecular modelling, and, more recently, drug discovery [41].
In the past decade, there has been immense progress in the development of algorithms and technology using mutations, which has transformed the field of protein design and engineering to attain tailormade proteins suitable for pharmaceutical and biotechnological applications. Among various in silico tools available, the Cologne University Protein Stability Analysis Tool (CUPSAT) [42], Site Directed Mutator (SDM) [43], PopMusic 2.1 [44], SNPeffect 4.0 [45], PolyPhen-2 [46], DUET [47], MAESTRO Web [48], DynaMut [49], and mCSM PPI2 [50] for mutational studies have been successfully used to evaluate stability change (stabilizing or destabilizing) and, after mutations, to predict the phenotypic consequence of missense variants. These tools are structure based, sequence based, and energy based, and combined features (statistical approach and/or machine learning methods, such as neural networks and support vector machines (SVM). These methods are fast, user-friendly, and reliable and promise to be invaluable in the development of proteins with a wide range of impactful applications. The applications of these tools extend from understanding the origin of diseases caused by misregulation of protein maintenance [51,52] to discriminating disease-associated mutations from non-disease mutations, studying drug-resistant mutations [53,54,55], and providing important structural and functional insights into designing new proteins [56,57,58,59]. To design DARPins as ERK2 inhibitors, a multi-point mutation approach, MAESTRO, was applied to the wild-type DARPin protein to identify the stabilizing mutation points, followed by validation of the binding energy of mutants employing MD simulations using MM-PBSA/GBSA protocols. The details of the results are discussed in the next section.

2. Results

2.1. Design and Prediction of New Inhibitors

The effect of mutations on the thermodynamic stability of DARPins (E40) was analysed using MAESTRO and the predictions were further analysed using other algorithms, and a comparison of ΔΔG predictions is shown in Table S1. MAESTRO is an easy and standalone program that provides different kinds of mutation experiments on single chains and protein complexes. The predictive power of this method is suggested to be reliable as it combines multiple linear regression (MLR), a neural network approach (NN), and a support vector machine (SVM) that allows to include additional information such as protein size or solvent accessibility. The mutation sensitivity profile of E40 is shown in Figure 2. Mutation points were selected based on MAESTRO suggestions, from which only randomized residues of DARPin E40 were chosen to undergo further investigation; refer to Figure 1B.

2.1.1. Evaluation of Specific Mutations

In total, seven single-point favourable mutations (S380, I389, D421, N422, A443, D454, and R455) were suggested by MAESTRO based on predicted change in stability and confidence estimation calculated in terms of ∆∆G and Cpred, respectively. Wild-type amino acids (AA) in DARPin E40 were mutated with 16 other available amino acids (except cysteine, proline, and glycine). Table 1 shows the suggested mutations with their respective ∆∆Gpred and Cpred.

2.1.2. Evaluation of Selected Mutants by MDS

In the light of these results, out of all the suggested mutation points from MAESTRO, a total of 13 mutation points meeting the criterion, i.e., ∆∆G < 0.0 (stabilizing) and Cpred values ~1 (highly reliable), were selected for further analysis. Their structures were modelled using the “Build mutant” protocol in DS Modeller [60], further optimized, and finally subjected to large-scale MD simulations to investigate the structural consequences of mutating residues. To study the forces interplay, trajectories of mutant complexes were analysed for binding free energy, per-residue binding free energy of complexes, and pairwise binding free energy of residues within 4 Å.
The trajectories obtained from the production run of 100 ns MD simulations of mutants were analysed for their binding free energy. Binding free energies of mutant complexes (A443D/ERK2, S380L/ERK2, A443N/ERK2, and N422T/ERK2) suggested better binding than E40/ERK2 (Table 2. In the next section, the decomposition of binding energy and important interacting residues with their H-bonds are discussed in detail. MAESTRO gave the highest match with other prediction techniques and has been used as a guide for the new design.
For ranking MM-GBSA, binding energy criteria are used and the energy contributions are Ggas and Gsolv. The polar and non-polar contributions are EGB and ESURF, respectively, for the GB calculations shown in Table 3. Although the electrostatic contribution from all mutants is big, it is mostly compensated by a large positive polar contribution (EGB), making the total polar contribution (EEL + EGB) positive and hence it is mostly the van der Waals term that contributes towards the total binding free energy.
To understand the residue contribution from DARPin and ERK2, a comparison was made for the decomposed free energy for each residue. Important binding residues are shown with their respective energies for all the three mutants, and E40/ERK2.DARPin residues start with “L,” while ERK2 residues start with “R” followed by a three-letter code of amino acids. From Figure 3 it is suggested that most of the important interactions from ERK2 come from the activation loop, specifically ARG180 contributing the most to the αG and L14 regions, while for all the DARPins (N3C), A443D, S380L, and D421W interactions come from all repeats (2,3), including L-ASP409, TRP413, and ASP421, but residues from C-cap (L-TYR444, ASP454, and PHE477) terminal contribute the most. Among all the complexes, A443D/ERK2 and S380L/ERK2 show stronger interactions with decomposed binding free energy <−5 kcal/mol for the above-mentioned residues.
To evaluate the effect of stabilizing mutations on E40 DARPin, the RMSD values of the position differences of backbone atoms between mutant and wild-type structures were calculated throughout the MD simulations. A comparison of the RMSD of E40/ERK2 with the other four complexes pE59/ERK2, A443D/ERK2, S380L/ERK2, and D421W/ERK2 is shown in Figure S1. According to the figures, all the simulations are well converged. According to MAESTRO suggestion, the mutations A443D, S380L, and D421W are stabilizing. After performing 100 ns simulation on the complexes A443D/ERK2, S380L/ERK2, and D421W/ERK2, the average RMSD values are ~3.19, 2.61, and 2.80 Å, respectively.

2.2. Exploring the Binding Mechanism of DARPins with ERK2

The main contribution towards the receptor-ligand affinity comes from non-covalent interactions. In the complex A443D/ERK2, alanine at position 443 on the DARPin loop is mutated to aspartate. Alanine contains a non-polar aliphatic R (CH3) group, which is small, being mutated to a polar aspartate that has a negatively charged or acidic R (CH2COO) group and a long side chain that can enhance electrostatic interactions, which can be seen from the electrostatic energy of −626.46 kcal/mol (Table 3). The second-best selected mutant of E40 that came out after MDS is S380L with a −7 kcal/mol dip in binding free energy (−56.74 kcal/mol) compared to the previous one (Table 2). Serine contains a hydroxymethyl group and is classified as a polar amino acid mutated to a positively charged and non-polar leucine with a side chain containing an isobutyl group. This mutation leads to a more negative electrostatic interaction energy of −673.80 kcal/mol than E40/ERK2 (Table 3).
To design a high binding inhibitor, it is required to understand interactions that differentiate it from the low binding inhibitors. For this purpose, the third mutant, D421W, having ∆∆G = −41.08 kcal/mol (Table 2), was selected. An aspartic acid at position 421 was mutated to tryptophan that contains a side chain indole, making it non-polar through its aromatic amino acid. Unlike the previous mutations that were located on the DARPin loops, this point mutation on the DARPin repeat has a negative effect on the binding affinity with a binding free energy of −41.08 kcal/mol compared to −49.50 kcal/mol for E40/ERK2.
Figure 4 displays pairwise decomposition free energy for residues within 4 Å. The DARPin interacts with the receptor mainly through the activation loop, alpha G and MAPK insertion regions [61]. The receptor-ligand interactions are mainly through hydrogen bonds (details of the type of hydrogen bonds are presented in Table 4). For all the three mutants, the hotspot residues show interactions mostly through cation-Π and salt bridges. The two most remarkable salt bridge interactions occur between ARG180-ASP454 and LYS220-ASP409 (Figure 5). In the mutant A443D, salt bridges ARG180-HH12–ASP454-OD1 and LYS220-HZ1–ASP409-OD2 possess quite a low polar interaction energy—(−8.76 and −8.25 kcal/mol) and (−8.98 and −8.27 kcal/mol), respectively; hence the total binding free energy of these residue pairs. There is a remarkable decrease in the energy of LYS220-HZ1–ASP409-OD2 here compared to E40/ERK2. For the next mutant complex S380L/ERK2, the pairwise interaction of 4 Å residues is almost like that in the previous case A443D/ERK2. Here, the interactions between receptor and ligand are also through hydrogen bonds (Table 5) but differ in their energies, with two strong salt bridge interactions between ARG180-HH22–ASP454-OD1 (−8.30 kcal/mol) and LYS220 HZ2-ASP409 OD2. A noticeable drop is observed in the decomposition energy of the receptor-ligand pair LYS220–ASP409 (−11.45 kcal/mol). For the third mutant complex, D421W/ERK2 also follows a pattern similar to that of the pairwise interaction of 4 Å residues like that in the previous cases of A443D and S380L. As in the previous cases, the interactions occur between receptor and ligand through hydrogen bonds, but are not as strong as the previous mutations. Except for a strong salt bridge interaction between ARG180-HH22–ASP454-OD1 (−8.30 kcal/mol), all other interactions have lower ∆∆G values compared to the other two mutants, while salt bridge LYS220 HZ2-ASP409 OD2 shows a remarkable increase (positive) in ∆∆G values. The main polar receptor-ligand interactions of the mutant complexes are shown in Figure 6.
Electrostatic interactions are a pivotal player in protein interactions and helpful in understanding intermolecular protein-protein interactions as they are long-range and have influence on charge molecules. Electrostatic potential maps are also called electrostatic potential energy maps or molecular electrical potential surfaces. These maps demonstrate the charge distributions of molecules in three dimensions and allow us to visualize the variably charged regions of a molecule. Knowledge about the charge distributions can be useful in determining how molecules interact with one another. Moreover, electrostatic forces help in fast recognizing the right partner among hundreds of thousands of candidates present in the intracellular environment of protein-protein complexes [62].
The electrostatic potential surfaces of the complexes E40/ERK2, A443D/ERK2, S80L/ERK2, and D421W/ERK2 were generated using the Adaptive Poisson-Boltzmann Solver (APBS) plugin in PyMOL [63] keeping the same orientation.
To study in detail the specific electrostatic potential of the interacting residues, the interface region in the map was zoomed and the potential was matched with the van der Waals and electrostatic energies of the interacting residues in receptor-ligand pairs.
In Figure 7, it can be observed that the interacting residues from receptors and ligands have opposite potentials, i.e., electropositive (blue) and electronegative (red). The electrostatic and van der Waals energy shown in the bottom left corner also corroborates the same. We know that van der Waals are weak forces or temporary attractions between electron-rich regions of one molecule and electron-poor regions of another. For the arginine residue from ARG180, a guanidino group is protonated to give the guanidinium form (-C-(NH2)2+), making arginine a charged, aliphatic amino acid; the lysine from LYS220 is also protonated −NH3+ at physiological pH, and the tyrosine from TYR222 contains 4-hydroxyphenylalanine that is neutral. These residue charges extend from positive to neutral, which is shown by blue to white regions in the maps. In order to interact, they need to recognize oppositely (negatively) charged residues (aspartate form, −COO−) in their vicinity that are ASP409, ASP454, and ASP475, shown by red regions in the maps. ARG180-ASP454 shows the strongest electrostatic energy because of the salt bridge formed between these residues. Although there exist various interactions between receptors and ligands, only important residues are shown in Figure 6.
In general, single-point mutation in ankyrin does not affect the secondary structure. We have predicted the protein secondary structure for the wild-type and the mutated ankyrin using PredictProtein (ROSTLAB, Technische Universität München, https://predictprotein.org/) [64] and SCRATCH protein predictor (Institute for Genomics and Bioinformatics, University of California, Irvine, USA, http://scratch.proteomics.ics.uci.edu/) [65]. Results are shown in Supplementary Material, Figure S2. It was observed that the secondary structures of mutant DARPins were unaffected.

3. Discussion

To design DARPins with higher binding affinity towards ERK2, the MAESTRO method was used to predict stability upon point mutations in N3C DARPin (E40). From the suggested mutation points, 13 stabilizing mutants were subjected to 100 ns simulations and then compared with the wild-type complex (E40/ERK2) in terms of binding free energy. After evaluation, two DARPins, A443D and S380L (showing higher binding affinity than E40), and D421W (showing lower binding than E40) were selected to further analyse their binding mechanism. Out of the suggested mutants, the key elements of the ankyrin mutants with high and low affinity towards ERK2 were compared in terms of binding free energy, decomposition free energy, and strength and type of hydrogen bonds formed between receptor and ligands. It was found that all the DARPins show interactions, mostly coming from their third repeat and C-cap terminal residues. The most persistent interactions in each structure were studied via a H-bond analysis, wherein it was found that the high binding affinity of all the N3C DARPins is mainly attributed to salt bridges (ARG180- ASP454 and LYS220-ASP409) along with various other DARPin–ERK2 interactions. Analysis of the energy decomposition shows that electrostatic and van der Waals interactions were the main contributors towards the total binding free energy of these systems. The binding free energy of the newly designed DARPins was improved up to ~20%. The binding affinity of these systems shows a trend of E40 < D421W < S380L < A443D.
The wild-type DARPin (E40) along with mutants interact (A443D, S380L, and D421W) through hydrogen bonds, mainly salt bridges and cation-Π. Moreover, the strength of the salt bridges formed between ARG180- ASP454 and LYS220-ASP409 plays a significant role in deciding the binding affinity of DARPins towards ERK2. For a DARPin to have a good binding affinity, both interactions must be strong, which can be shown by the order of decomposition free energy of these receptor-ligand pairs and the total binding free energy of the complex. Electrostatic potential surface analysis revealed that mutants that can generate an electronegative potential near the binding interface, show good binding with ERK2. It was also observed that mutations on the loops of DARPins extend better binding compared to that on a DARPin repeat. All this information can be used in the design of new DARPin inhibitors against ERK2.

4. Materials and Methods

The starting structure for the present study is the X-ray crystal structure of DARPin E40 complexed with ERK2 (PDB ID: 3ZU7) [66], designated as E40/ERK2 for further reference in this study. For the complex (E40/ERK2), chain A of ERK2 and chain C of DARPin (E40) were taken as the initial structure for the MD simulations (Figure 8).
MD simulations were performed using the GPU version of Particle Mesh Ewald Molecular Dynamics (PMEMD.CUDA) from AMBER14 [67]. At the molecular level, physical forces were implemented using ff14SB [68] protein force field to carry out MD simulations. The structures were solvated (using the tLeap module implemented in AMBER). The water model TIP3P [69] was used, wherein a cubic box of water extends at least 10 Å from the solute in each direction. A cut-off distance of 15 Å was used to compute the non-bonded interactions (electrostatic interactions and van der Waals interactions). For each of the complexes, a 100 ns long simulation was performed. To minimize the edge effects, all simulations were performed under periodic boundary conditions, and to treat long-range electrostatics, the particle mesh Ewald method was used [70]. To relax the system prior to MD simulation, the complexes were minimized using a series of steepest descent (SD) and conjugated gradient (CG) under the sander module of the AMBER14 program. During the simulation, the system was heated gradually over a period of 60 ps from 0 to 310 K (biological temperature), and a force constant of 5 kcal/mol Å2 was applied to restrain the atomic position. A constant pressure of 1 atm for 200 ps (NPT) was applied under Langevin dynamics, followed by a 40 ps volume-constant period (NVT) at a force constant of 2.5 kcal/mol Å2, which was maintained and followed by 100 ps dynamics at a force constant of 1.25 kcal/mol Å2. Finally, unrestrained production runs were performed for 100 ns, wherein no force was applied on any protein atoms in the NVT ensemble at a constant temperature of 310 K (biological temperature). For all analyses, 500 snapshots were taken from the last 5 ns of the simulation (96–100 ns). To check the system equilibrium, root-mean-square deviation (RMSD) of all backbone atoms was performed using the cpptraj module [71] incorporated in AmberTools 15 and compared to the starting structure. All simulations were carried out under periodic boundary conditions [72]. To treat long-range electrostatics, the particle mesh Ewald method was used [73,74,75], and the SHAKE algorithm was employed to constrain bond lengths involving hydrogen atoms. A 2 fs time step was set up while the trajectory was recorded every 0.1 ps. To relax the system prior to MD simulation, a series of steepest descent (SD) along with conjugated gradient (CG) minimizations with a total of 500 steps each was performed. Post-processing of the trajectories was performed using the MM-PBSA/GBSA protocol [76].

Supplementary Materials

The following are available online: Figure S1: Comparison of the RMSD values of mutants and E40/ERK2; Figure S2: Comparison of secondary structures of DARPin E40 and mutants (A443D, S380L, and D421W) obtained by using (A) PredictProtein and (B) SCRATCH protein predictor web servers. (For PredictProtein, helixes are shown in blue and others in yellow; for the SCRATCH protein predictor, H is for helix and C is for others.); Table S1: Comparison of predicted ∆∆G (kcal/mol) mutants of 3ZU7 from different web-based algorithms.

Author Contributions

Conceptualization, V.G. and V.S.L.; funding acquisition, P.N.; investigation, V.G.; methodology, V.G.; resources, S.M.Z. and V.S.L.; supervision, S.M.Z., N.A.R., and V.S.L.; writing—original draft, V.G. All authors have read and agreed to the published version of the manuscript.

Funding

This research is conducted with the support by Ministry of Higher Education (MOHE) through Fundamental Research Grant Scheme (FRGS/1/2019/STG03/UM/02/8), computer facility by Data Intensive Computing Centre (DICC) Universiti Malaya, the Distinguished Research Professor Grant (NRCT 808/2563) of the National Research Council of Thailand, the Office of National Higher Education Science Research and Innovation Policy Council, the Program Management Unit for Human Resources & Institutional Development, Research and Innovation [B05F630102], the National Science and Technology Development Agency (NSTDA), the National Research University project under Thailand’s Office of the Higher Education Commission, the Center of Excellence for Innovation in Analytical Science and Technology (I-ANALY-S-T), and the Visiting Professor Grant from Chiang Mai University, Thailand (R000026614).

Data Availability Statement

The data presented in this study are available within the article or supplementary materials.

Conflicts of Interest

The authors declare no conflict of interest.

Sample Availability

Samples of the compounds are available upon request.

References

  1. Roskoski, R. ERK1/2 MAP kinases: Structure, function, and regulation. Pharmacol. Res. 2012, 66, 105–143. [Google Scholar] [CrossRef]
  2. Pearson, G.; Robinson, F.; Gibson, T.B.; Xu, B.E.; Karandikar, M.; Berman, K.; Cobb, M.H. Mitogen-activated protein (MAP) kinase pathways: Regulation and physiological functions. Endocr. Rev. 2001, 22, 153–183. [Google Scholar]
  3. Bhullar, K.S.; Lagarón, N.O.; McGowan, E.M.; Parmar, I.; Jha, A.; Hubbard, B.P.; Rupasinghe, H.P.V. Kinase-targeted cancer therapies: Progress, challenges and future directions. Mol. Cancer 2018, 17, 1–20. [Google Scholar] [CrossRef]
  4. Dreier, B.; Plückthun, A. Rapid Selection of High-Affinity Antibody scFv Fragments Using Ribosome Display. Methods Mol. Biol. 2018, 1827, 235–268. [Google Scholar] [PubMed]
  5. Steiner, D.; Forrer, P.; Plückthun, A. Efficient selection of DARPins with sub-nanomolar affinities using SRP phage display. J. Mol. Biol. 2008, 382, 1211–1227. [Google Scholar] [CrossRef]
  6. Hartmann, J.; Münch, R.C.; Freiling, R.T.; Schneider, I.C.; Dreier, B.; Samukange, W.; Koch, J.; Seeger, M.A.; Plückthun, A.; Buchholz, C.J. A Library-Based Screening Strategy for the Identification of DARPins as Ligands for Receptor-Targeted AAV and Lentiviral Vectors. Mol. Ther. Methods Clin. Dev. 2018, 10, 128–143. [Google Scholar] [CrossRef] [PubMed]
  7. Smithwick, E.; Stewart, M.W. Designed Ankyrin Repeat Proteins: A Look at their Evolving Use in Medicine with a Focus on the Treatment of Chorioretinal Vascular Disorders. Antiinflamm. Antiallergy Agents Med. Chem. 2017, 16, 33–45. [Google Scholar] [CrossRef]
  8. Proshkina, G.; Deyev, S.; Ryabova, A.; Tavanti, F.; Menziani, M.C.; Cohen, R.; Katrivas, L.; Kotlyar, A. DARPin_9-29-Targeted Mini Gold Nanorods Specifically Eliminate HER2-Overexpressing Cancer Cells. ACS Appl. Mater Interfaces 2019, 11, 34645–34651. [Google Scholar] [CrossRef] [PubMed]
  9. Shilova, O.N.; Deyev, S.M. DARPins: Promising Scaffolds for Theranostics. Acta. Nat. 2019, 11, 42–53. [Google Scholar] [CrossRef]
  10. Stumpp, M.T.; Dawson, K.M.; Binz, H.K. Beyond Antibodies: The DARPin® Drug Platform. BioDrugs 2020, 34, 423–433. [Google Scholar] [CrossRef]
  11. Boersma, Y.L.; Plückthun, A. DARPins and other repeat protein scaffolds: Advances in engineering and applications. Curr. Opin. Biotechnol. 2011, 22, 849–857. [Google Scholar] [CrossRef]
  12. Forrer, P.; Binz, H.K.; Stumpp, M.T.; Plückthun, A. Consensus Design of Repeat Proteins. ChemBioChem 2004, 5, 183–189. [Google Scholar] [CrossRef]
  13. Forrer, P.; Stumpp, M.T.; Binz, H.K.; Plückthun, A. A novel strategy to design binding molecules harnessing the modular nature of repeat proteins. FEBS Lett. 2003, 39, 2–6. [Google Scholar] [CrossRef]
  14. Kohl, A.; Binz, H.K.; Forrer, P.; Stumpp, M.T.; Plückthun, A.; Grütter, M.G. Designed to be stable: Crystal structure of a consensus ankyrin repeat protein. Proc. Natl. Acad. Sci. USA 2003, 100, 1700–1705. [Google Scholar] [CrossRef] [PubMed]
  15. Gorina, S.; Pavletich, N.P. Structure of the p53 tumor suppressor bound to the ankyrin and SH3 domains of 53BP2. Science 1996, 274, 1001–1005. [Google Scholar] [CrossRef]
  16. Sedgwick, S.G.; Smerdon, S.J. The ankyrin repeat: A diversity of interactions on a common structural framework. Trends Biochem. Sci. 1999, 24, 311–316. [Google Scholar] [CrossRef]
  17. Binz, H.K.; Stumpp, M.T.; Forrer, P.; Amstutz, P.; Plückthun, A. Designing repeat proteins: Well-expressed, soluble and stable proteins from combinatorial libraries of consensus ankyrin repeat proteins. J. Mol. Biol. 2003, 332, 489–503. [Google Scholar] [CrossRef]
  18. Leader, B.; Baca, Q.J.; Golan, D.E. Protein therapeutics: A summary and pharmacological classification. Nat. Rev. Drug Discov. 2008, 7, 21–39. [Google Scholar] [CrossRef]
  19. Gautieri, A.; Beeg, M.; Gobbi, M.; Rigoldi, F.; Colombo, L.; Salmona, M. The Anti-Amyloidogenic Action of Doxycycline: A Molecular Dynamics Study on the Interaction with Aβ42. Int. J. Mol. Sci. 2019, 20, 4641. [Google Scholar] [CrossRef] [PubMed]
  20. Song, X.; Jensen, M.Ø.; Jogini, V.; Stein, R.A.; Lee, C.H.; Mchaourab, H.S.; Shaw, D.E.; Gouaux, E. Mechanism of NMDA receptor channel block by MK-801 and memantine. Nature 2018, 556, 515–519. [Google Scholar] [CrossRef] [PubMed]
  21. Sun, N.; Yuan, C.; Ma, X.; Wang, Y.; Gu, X.; Fu, W. Molecular Mechanism of Action of RORγt Agonists and Inverse Agonists: Insights from Molecular Dynamics Simulation. Molecules 2018, 23, 3181. [Google Scholar] [CrossRef]
  22. Wang, H.; Gao, Z.; Song, P.; Hu, B.; Wang, J.; Cheng, M. Molecular dynamics simulation and QM/MM calculation reveal the selectivity mechanism of type I 1/2 kinase inhibitors: The effect of intramolecular H-bonds and conformational restriction for improved selectivity. Phys. Chem. Chem. Phys. 2019, 21, 24147–24164. [Google Scholar] [CrossRef]
  23. Vila-Viçosa, D.; Victor, B.L.; Ramos, J.; Machado, D.; Viveiros, M.; Switala, J.; Loewen, P.C.; Leitão, R.; Martins, F.; Machuqueiro, M. Insights on the Mechanism of Action of INH-C10 as an Antitubercular Prodrug. Mol. Pharm. 2017, 14, 4597–4605. [Google Scholar] [CrossRef]
  24. Barros, E.P.; Schiffer, J.M.; Vorobieva, A.; Dou, J.; Baker, D.; Amaro, R.E. Improving the Efficiency of Ligand-Binding Protein Design with Molecular Dynamics Simulations. J. Chem. Theory Comput. 2019, 15, 5703–5715. [Google Scholar] [CrossRef] [PubMed]
  25. Chew, T.A.; Orlando, B.J.; Zhang, J.; Latorraca, N.R.; Wang, A.; Hollingsworth, S.A.; Chen, D.H.; Dror, R.O.; Liao, M.; Feng, L. Structure and mechanism of the cation-chloride cotransporter NKCC1. Nature 2019, 572, 488–492. [Google Scholar] [CrossRef]
  26. Ferreira, J.V.; Capello, T.M.; Siqueira, L.J.; Lago, J.H.; Caseli, L. Mechanism of Action of Thymol on Cell Membranes Investigated through Lipid Langmuir Monolayers at the Air-Water Interface and Molecular Simulation. Langmuir 2016, 32, 3234–3241. [Google Scholar] [CrossRef] [PubMed]
  27. Yang, T.; Zhang, W.; Cheng, J.; Nie, Y.; Xin, Q.; Yuan, S.; Dou, Y. Formation Mechanism of Ion Channel in Channelrhodopsin-2: Molecular Dynamics Simulation and Steering Molecular Dynamics Simulations. Int. J. Mol. Sci. 2019, 20, 3780. [Google Scholar] [CrossRef] [PubMed]
  28. Behmard, E.; Najafi, A.; Ahmadi, A. Understanding the resistance mechanism of penicillin binding protein 1a mutant against cefotaxime using molecular dynamic simulation. J. Biomol. Struct. Dyn. 2019, 37, 741–749. [Google Scholar] [CrossRef]
  29. Chen, Q.; Cheng, X.; Wei, D.; Xu, Q. Molecular dynamics simulation studies of the wild type and E92Q/N155H mutant of Elvitegravir-resistance HIV-1 integrase. Interdiscip. Sci. 2015, 7, 36–42. [Google Scholar] [PubMed]
  30. Cloete, R.; Kapp, E.; Joubert, J.; Christoffels, A.; Malan, S.F. Molecular modelling, and simulation studies of the Mycobacterium tuberculosis multidrug efflux pump protein Rv1258c. PLoS ONE 2018, e0207605. [Google Scholar] [CrossRef]
  31. Ge, Y.; Wu, J.; Xia, Y.; Yang, M.; Xiao, J.; Yu, J. Molecular dynamics simulation of the complex PBP-2x with drug cefuroxime to explore the drug resistance mechanism of Streptococcus suis R61. PLoS ONE 2012, 7, e35941. [Google Scholar] [CrossRef] [PubMed]
  32. Ul Haq, F.; Abro, A.; Raza, S.; Liedl, K.R.; Azam, S.S. Molecular dynamics simulation studies of novel β-lactamase inhibitor. J. Mol. Graph. Model. 2017, 74, 143–152. [Google Scholar] [CrossRef] [PubMed]
  33. Aci-Sèche, S.; Ziada, S.; Braka, A.; Arora, R.; Bonnet, P. Advanced molecular dynamics simulation methods for kinase drug discovery. Future Med. Chem. 2016, 8, 545–566. [Google Scholar] [CrossRef] [PubMed]
  34. Ahuja, L.G.; Aoto, P.C.; Kornev, A.P.; Veglia, G.; Taylor, S.S. Dynamic allostery-based molecular workings of kinase:peptide complexes. Proc. Natl. Acad. Sci. USA 2019, 116, 15052–15061. [Google Scholar] [CrossRef]
  35. Bao, Y.; Zhou, L.; Dai, D.; Zhu, X.; Hu, Y.; Qiu, Y. Discover potential inhibitors for PFKFB3 using 3D-QSAR, virtual screening, molecular docking and molecular dynamics simulation. J. Recept. Signal Transduct. Res. 2018, 38, 413–431. [Google Scholar] [CrossRef]
  36. Li, K.; Zhu, J.; Xu, L.; Jin, J. Rational Design of Novel Phosphoinositide 3-Kinase Gamma (PI3Kγ) Selective Inhibitors: A Computational Investigation Integrating 3D-QSAR, Molecular Docking and Molecular Dynamics Simulation. Chem. Biodivers. 2019, 16, e1900105. [Google Scholar] [CrossRef]
  37. Liu, J.; Lu, Y.; Li, G.; Xiao, M.; Yang, G.; Pan, Y. Elucidation the binding mechanism of Nelumbo nucifera-derived isoquinoline alkaloids as Rho-kinase 1 inhibitors by molecular docking and dynamic simulation. J. Biomol. Struct Dyn. 2021, 39, 379–394. [Google Scholar] [CrossRef]
  38. Zou, Y.; Ewalt, J.; Ng, H.L. Recent Insights from Molecular Dynamics Simulations for G Protein-Coupled Receptor Drug Discovery. Int. J. Mol. Sci. 2019, 20, 4237. [Google Scholar] [CrossRef]
  39. Ma, C.; Hu, B.; Zhang, L.; Zhao, Y.; Wang, M.; Wang, J.; Cheng, M. Computational investigation of the antagonism effect towards GluN2B-Containing NMDA receptor: Combined ligand-based and target-based approach. J. Mol. Graph. Model. 2019, 86, 95–105. [Google Scholar] [CrossRef]
  40. Michael, E.; Polydorides, S.; Simonson, T.; Archontis, G. Hybrid MC/MD for protein design. J. Chem. Phys. 2020, 153, 054113. [Google Scholar] [CrossRef]
  41. De Vivo, M.; Masetti, M.; Bottegoni, G.; Cavalli, A. Role of Molecular Dynamics and Related Methods in Drug Discovery. J. Med. Chem. 2016, 59, 4035–4061. [Google Scholar] [CrossRef]
  42. Parthiban, V.; Gromiha, M.M.; Schomburg, D. CUPSAT: Prediction of protein stability upon point mutations. Nucleic Acids Res. 2006, 1, W239–W242. [Google Scholar] [CrossRef] [PubMed]
  43. Worth, C.L.; Preissner, R.; Blundell, T.L. SDM—A server for predicting effects of mutations on protein stability and malfunction. Nucleic Acids Res. 2011, 39, W215–W222. [Google Scholar] [CrossRef] [PubMed]
  44. Dehouck, Y.; Kwasigroch, J.M.; Gilis, D.; Rooman, M. PoPMuSiC 2.1: A web server for the estimation of protein stability changes upon mutation and sequence optimality. BMC Bioinformat. 2011, 13, 151. [Google Scholar] [CrossRef]
  45. De Baets, G.; Van Durme, J.; Reumers, J.; Maurer-Stroh, S.; Vanhee, P.; Dopazo, J.; Schymkowitz, J.; Rousseau, F. SNPeffect 4.0: On-line prediction of molecular and structural effects of protein-coding variants. Nucleic Acids Res. 2012, 40, D935–D939. [Google Scholar] [CrossRef] [PubMed]
  46. Adzhubei, I.; Jordan, D.M.; Sunyaev, S.R. Predicting functional effect of human missense mutations using PolyPhen-2. Curr. Protoc. Hum. Genet. 2013, 76, 7–20. [Google Scholar] [CrossRef] [PubMed]
  47. Pires, D.E.; Ascher, D.B.; Blundell, T.L. mCSM: Predicting the effects of mutations in proteins using graph-based signatures. Bioinformatics 2014, 30, 335–342. [Google Scholar] [CrossRef] [PubMed]
  48. Laimer, J.; Hofer, H.; Fritz, M.; Wegenkittl, S.; Lackner, P. MAESTRO-multi agent stability prediction upon point mutations. BMC Bioinformat. 2015, 16, 116. [Google Scholar] [CrossRef] [PubMed]
  49. Rodrigues, C.H.; Pires, D.E.; Ascher, D.B. DynaMut: Predicting the impact of mutations on protein conformation, flexibility and stability. Nucleic Acids Res. 2018, 2, W350–W355. [Google Scholar] [CrossRef]
  50. Rodrigues, C.H.M.; Myung, Y.; Pires, D.E.V.; Ascher, D.B. mCSM-PPI2: Predicting the effects of mutations on protein-protein interactions. Nucleic Acids Res. 2019, 2, W338–W344. [Google Scholar] [CrossRef]
  51. Birolo, G.; Benevenuta, S.; Fariselli, P.; Capriotti, E.; Giorgio, E.; Sanavia, T. Protein Stability Perturbation Contributes to the Loss of Function in Haploinsufficient Genes. Front. Mol. Biosci. 2021, 1, 620793. [Google Scholar] [CrossRef]
  52. Gerasimavicius, L.; Liu, X.; Marsh, J.A. Identification of pathogenic missense mutations using protein stability predictors. Sci. Rep. 2020, 21, 15387. [Google Scholar] [CrossRef] [PubMed]
  53. Powers, E.T.; Morimoto, R.I.; Dillin, A.; Kelly, J.W.; Balch, W.E. Biological and chemical approaches to diseases of proteostasis deficiency. Annu. Rev. Biochem. 2009, 78, 959–991. [Google Scholar] [CrossRef]
  54. Karmakar, M.; Rodrigues, C.H.M.; Holt, K.E.; Dunstan, S.J.; Denholm, J.; Ascher, D.B. Empirical ways to identify novel Bedaquiline resistance mutations in AtpE. PLoS ONE 2019, 29, 14. [Google Scholar] [CrossRef]
  55. Pines, G.; Fankhauser, R.G.; Eckert, C.A. Predicting Drug Resistance Using Deep Mutational Scanning. Molecules 2020, 25, 2265. [Google Scholar] [CrossRef]
  56. Portelli, S.; Phelan, J.E.; Ascher, D.B.; Clark, T.G.; Furnham, N. Understanding molecular consequences of putative drug resistant mutations in Mycobacterium tuberculosis. Sci. Rep. 2018, 8, 15356. [Google Scholar] [CrossRef]
  57. Fischer, A.; Seitz, T.; Lochner, A.; Sterner, R.; Merkl, R.; Bocola, M. A fast and precise approach for computational saturation mutagenesis and its experimental validation by using an artificial (βα)8-barrel protein. Chembiochem 2011, 12, 1544–1550. [Google Scholar] [CrossRef] [PubMed]
  58. Munir, A.; Vedithi, S.C.; Chaplin, A.K.; Blundell, T.L. Genomics, Computational Biology and Drug Discovery for Mycobacterial Infections: Fighting the Emergence of Resistance. Front. Genet. 2020, 11, 965. [Google Scholar] [CrossRef] [PubMed]
  59. Xin, L.; Yu, H.; Hong, Q.; Bi, X.; Zhang, X.; Zhang, Z.; Kong, Z.; Zheng, Q.; Gu, Y.; Zhao, Q.; et al. Identification of Strategic Residues at the Interface of Antigen-Antibody Interactions by In Silico Mutagenesis. Interdiscip. Sci. 2018, 10, 438–448. [Google Scholar] [CrossRef]
  60. BIOVIA, D.S. Dassault Systèmes BIOVIA, (Version 9.24); Dassault Systèmes: San Diego, CA, USA, 2016. [Google Scholar]
  61. Gautam, V.; Chong, W.L.; Chin, S.P.; Zain, S.M.; Rahman, N.A.; Vao-soongnern, V.; Lee, V.S. Loop dynamics behind the affinity of DARPins towards ERK2: Molecular dynamics simulations (MDs) and elastic network model (ENM). J. Mol. Liq. 2019, 274, 612–620. [Google Scholar] [CrossRef]
  62. Vascon, F.; Gasparotto, M.; Giacomello, M.; Cendron, L.; Bergantino, E.; Filippini, F.; Righetto, I. Protein electrostatics: From computational and structural analysis to discovery of functional fingerprints and biotechnological design. Comput. Struct. Biotechnol. J. 2020, 18, 1774–1789. [Google Scholar] [CrossRef]
  63. The PyMOL Molecular Graphics System, Version 2.0 Schrödinger, LLC.
  64. Yachdav, G.; Kloppmann, E.; Kajan, L.; Hecht, M.; Goldberg, T.; Hamp, T.; Hönigschmid, P.; Schafferhans, A.; Roos, M.; Bernhofer, M.; et al. PredictProtein—an open resource for online prediction of protein structural and functional features. Nucleic Acids Res. 2014, 42, W337–W343. [Google Scholar] [CrossRef] [PubMed]
  65. Cheng, J.; Randall, A.Z.; Sweredoski, M.J.; Baldi, P. SCRATCH: A protein structure and structural feature prediction server. Nucleic Acids Res. 2005, 33, W72–W76. [Google Scholar] [CrossRef]
  66. Kummer, L.; Parizek, P.; Rube, P.; Millgramm, B.; Prinz, A.; Mittl, P.R.; Kaufholz, M.; Zimmermann, B.; Herberg, F.W.; Plückthun, A. Structural and functional analysis of phosphorylation-specific binders of the kinase ERK from designed ankyrin repeat protein libraries. Proc. Natl. Aacd. Sci. USA 2012, 109, E2248–E2257. [Google Scholar] [CrossRef]
  67. Case, D.A.; Babin, V.; Berryman, J.T.; Betz, R.M.; Cerutti, D.; Cheatham, T.E., III; Darden, T.A.; Duke, R.E.; Kollman, P.A. AMBER 14; University of California: San Francisco, CA, USA, 2014. [Google Scholar]
  68. Pérez, A.; Marchán, I.; Svozil, D.; Sponer, J.; Cheatham, T.E., III; Laughton, C.A.; Orozco, M. Refinement of the AMBER force field for nucleic acids: Improving the description of alpha/gamma conformers. Biophys. J. 2007, 92, 3817–3829. [Google Scholar] [CrossRef] [PubMed]
  69. MacKerell, A.D.; Bashford, D.; Bellott, M.; Dunbrack, R.L.; Evanseck, J.D.; Field, M.J.; Fischer, S.; Karplus, M. All-atom empirical potential for molecular modeling and dynamics studies of proteins. J. Phys. Chem. B. 1998, 102, 3586–3616. [Google Scholar] [CrossRef]
  70. Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N⋅log(N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089–10092. [Google Scholar] [CrossRef]
  71. Roe, D.R.; Cheatham, T.E., III. PTRAJ and CPPTRAJ: Software for Processing and Analysis of Molecular Dynamics Trajectory Data. J. Chem. Theory Comput. 2013, 9, 3084–3095. [Google Scholar] [CrossRef]
  72. Linse, B.; Linse, P. Tuning the smooth particle mesh Ewald sum: Application on ionic solutions and dipolar fluids. J. Chem. Phys. 2014, 141, 184114. [Google Scholar] [CrossRef]
  73. Darden, T.; Perera, L.; Li, L.; Pedersen, L. New tricks for modelers from the crystallography toolkit: The particle mesh Ewald algorithm and its use in nucleic acid simulations. Structure 1999, 15, R55–R60. [Google Scholar] [CrossRef]
  74. Kollman, P.A.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Lee, T.; Duan, Y.; Wang, W.; et al. Calculating structures and free energies of complex molecules: Combining molecular mechanics and continuum models. Acc. Chem. Res. 2000, 33, 889–897. [Google Scholar] [CrossRef] [PubMed]
  75. Kuhn, B.; Kollman, P.A. Binding of a diverse set of ligands to avidin and streptavidin: An accurate quantitative prediction of their relative affinities by a combination of molecular mechanics and continuum solvent models. J. Med. Chem. 2000, 43, 3786–3791. [Google Scholar] [CrossRef] [PubMed]
  76. Genheden, S.; Ryde, U. The MM/PBSA and MM/GBSA methods to estimate ligand-binding affinities. Expert Opin. Drug Discov. 2015, 10, 449–461. [Google Scholar] [CrossRef] [PubMed]
Molecules 26 04540 g001 550
Figure 1. (A) Architecture of DARPins (N3C) showing C-cap, N-cap, and internal repeats. (B) Consensus design of a DARPin repeat containing 33 amino acids showing framework (black) and randomized residues (red).
Figure 1. (A) Architecture of DARPins (N3C) showing C-cap, N-cap, and internal repeats. (B) Consensus design of a DARPin repeat containing 33 amino acids showing framework (black) and randomized residues (red).
Molecules 26 04540 g001
Molecules 26 04540 g002 550
Figure 2. Mutation-sensitive profile of the complex E40/ERK2 obtained from MAESTRO. Confidence estimation and prediction error on multi-point mutations are shown. Blue and red bars show stabilizing and destabilizing mutations, respectively.
Figure 2. Mutation-sensitive profile of the complex E40/ERK2 obtained from MAESTRO. Confidence estimation and prediction error on multi-point mutations are shown. Blue and red bars show stabilizing and destabilizing mutations, respectively.
Molecules 26 04540 g002
Molecules 26 04540 g003 550
Figure 3. Comparison between the per-residue free energy decomposition of E40/ERK2 and mutants. DARPin residues start with “L”, while ERK2 residues start with “R” followed by a three-letter code of amino acids.
Figure 3. Comparison between the per-residue free energy decomposition of E40/ERK2 and mutants. DARPin residues start with “L”, while ERK2 residues start with “R” followed by a three-letter code of amino acids.
Molecules 26 04540 g003
Molecules 26 04540 g004 550
Figure 4. Comparison of pairwise decomposition free energy of complexes. Clockwise: E40/ERK2, A443D/ERK2, D421W/ERK2, and S380L/ERK2. Residues of receptor-ligand pairs and decomposition free energies in kcal/mol are shown on the X and Y axes, respectively. Contributions of van der Waals, non-polar solvation, polar, and total energies are shown with different colour bars.
Figure 4. Comparison of pairwise decomposition free energy of complexes. Clockwise: E40/ERK2, A443D/ERK2, D421W/ERK2, and S380L/ERK2. Residues of receptor-ligand pairs and decomposition free energies in kcal/mol are shown on the X and Y axes, respectively. Contributions of van der Waals, non-polar solvation, polar, and total energies are shown with different colour bars.
Molecules 26 04540 g004
Molecules 26 04540 g005 550
Figure 5. The main salt bridges formed between DARPins and ERK2. (A) ARG180-NH2–ASP454-OD1 and (B) LYS220 NZ2-ASP409 OD2. Salt bridges are shown by black dots.
Figure 5. The main salt bridges formed between DARPins and ERK2. (A) ARG180-NH2–ASP454-OD1 and (B) LYS220 NZ2-ASP409 OD2. Salt bridges are shown by black dots.
Molecules 26 04540 g005
Molecules 26 04540 g006 550
Figure 6. Polar contacts of the complexes of ERK2; residues shown in yellow with (A) E40 in green, (B) S380L in blue, (C) A443D in orange, and (D) D421W in pink. Polar interactions are shown by black dotted lines.
Figure 6. Polar contacts of the complexes of ERK2; residues shown in yellow with (A) E40 in green, (B) S380L in blue, (C) A443D in orange, and (D) D421W in pink. Polar interactions are shown by black dotted lines.
Molecules 26 04540 g006
Molecules 26 04540 g007 550
Figure 7. Electrostatic potential surface of interacting residues (enlarged) in (A) A443D/ERK2 (B) S380L/ERK2, and (C) D421W/ERK2. Highly electronegative regions shown by red while electropositive regions are shown by blue color. Respective van der Waals and electrostatic energies for the receptor-ligand pairs are shown in bottom left corner. Red indicates highly electronegative, while blue shows electropositive residue surfaces (see bar). The areas of high electrostatic potential are found near the interface of DARPin and ERK2 in all the structures.
Figure 7. Electrostatic potential surface of interacting residues (enlarged) in (A) A443D/ERK2 (B) S380L/ERK2, and (C) D421W/ERK2. Highly electronegative regions shown by red while electropositive regions are shown by blue color. Respective van der Waals and electrostatic energies for the receptor-ligand pairs are shown in bottom left corner. Red indicates highly electronegative, while blue shows electropositive residue surfaces (see bar). The areas of high electrostatic potential are found near the interface of DARPin and ERK2 in all the structures.
Molecules 26 04540 g007
Molecules 26 04540 g008 550
Figure 8. Biological assembly of E40/ERK2 (PDB ID:3ZU7). Ribbon diagrams of DARPin and ERK2 are shown in orange and green, respectively.
Figure 8. Biological assembly of E40/ERK2 (PDB ID:3ZU7). Ribbon diagrams of DARPin and ERK2 are shown in orange and green, respectively.
Molecules 26 04540 g008
Table
Table 1. Specific mutation evaluation using the MAESTRO algorithm.
Table 1. Specific mutation evaluation using the MAESTRO algorithm.
MutantsΔΔGpredCpredMutantsΔΔGpredCpredMutantsΔΔGpredCpredMutantsΔΔGpredCpredMutantsΔΔGpredCpredMutantsΔΔGpredCpredMutantsΔΔGpredCpred
S380(A)0.0380.959N422(A)−0.5170.907A443(R)0.0670.807D454(A)−0.0350.933R455(A)−0.1460.922D421(A)−0.0310.926I389(Y)−0.7090.861
S380(R)0.4020.76N422(R)−0.050.851A443(N)−0.2320.92D454(R)0.1730.801R455(N)−0.270.915D421(A)0.1450.808I389(W)-0.6720.848
S380(N)−0.0070.946N422(D)−0.6970.835A443(D)−0.3250.947D454(N)−0.0940.906R455(D)−0.1710.934D421(A)−0.1150.909I389(V)−0.450.855
S380(D)0.0840.944N422(E)−0.7390.832A443(E)−0.3850.948D454(Q)−0.1250.909R455(E)−0.3050.931D421(A)−0.4610.882I389(T)−0.5220.846
S380(E)0.0650.937N422(Q)−0.4410.872A443(Q)−0.2780.925D454(H)0.0160.894R455(Q)−0.4390.895D421(A)−0.1280.913I389(S)−0.3780.84
S380(Q)−0.0750.915N422(H)−0.3590.898A443(H)0.1270.896D454(I)−0.2870.886R455(H)−0.5170.864D421(A)0.0230.901I389(R)−0.0590.881
S380(H)−0.1020.878N422(I)−0.4840.912A443(I)0.0560.905D454(L)−0.3520.867R455(I)−0.4330.903D421(A)−0.240.913I389(Q)−0.6440.858
S380(M)−0.2270.896N422(L)−0.4120.902A443(L)0.0410.902D454(K)0.4110.79R455(L)−0.4010.92D421(A)−0.2550.908I389(N)−0.5560.84
S380(F)−0.3460.898N422(K)−0.0670.909A443(K)0.0280.848D454(M)-0.3350.876R455(K)−0.5340.881D421(A)0.5030.795I389(M)−0.5890.871
S380(W)−0.4660.86N422(M)−0.4840.895A443(M)−0.1330.922D454(F)−0.3470.888R455(M)−0.4730.885D421(A)−0.2560.894I389(L)−0.4490.874
S380(T)−0.0020.944N422(F)−0.5120.879A443(F)−0.0520.914D454(S)0.0720.93R455(F)−0.5020.889D421(A)−0.3510.891I389(K)−0.0540.893
S380(Y)−0.3920.872N422(S)−0.5360.899A443(S)−0.2460.931D454(T)−0.1160.924R455(S)−0.1550.917D421(A)0.0330.941I389(H)−0.4340.875
S380(V)−0.2920.905N422(T)−0.6080.902A443(T)0.0390.909D454(W)−0.3170.873R455(T)−0.30.925D421(A)−0.4280.872I389(F)−0.7040.859
S380(K)0.2780.807N422(W)−0.7640.842A443(W)−0.0060.901D454(Y)−0.3520.888R455(W)−0.5680.879D421(A)−0.0210.955I389(E)-0.7930.83
S380(L)−0.2850.907N422(Y)−0.650.866A443(Y)−0.0010.905D454(V)−0.2150.921R455(Y)−0.5860.879D421(A)−0.4010.876I389(D)−0.8270.834
S380(I)−0.2090.904N422(V)−0.4790.899A443(V)0.0540.909D454(E)−0.3570.891R455(V)−0.3280.928D421(A)−0.1630.919I389(A)−0.2940.854
ΔΔGpred, total predicted change in stability (kcal/mol). ΔΔGpred < 0.0 indicates a stabilizing mutation. Cpred, confidence estimation, given as a value between 0.0 (not reliable) and 1.0 (highly reliable).
Table
Table 2. Comparison of the binding free energy (±SEM) of mutants obtained by the MM-PBSA/GBSA method.
Table 2. Comparison of the binding free energy (±SEM) of mutants obtained by the MM-PBSA/GBSA method.
MutantsMM-PBSA (kcal/mol)MM-GBSA (kcal/mol)
A443D−91.51 ± 0.42−59.86 ± 0.29
S380L−94.90 ± 0.34−56.74 ± 0.27
A443N−80.29 ± 0.34−51.59 ± 0.26
N422A−66.44 ± 0.39−38.87 ± 0.27
N422I−64.34 ± 0.34−37.49 ± 0.27
N422T−86.74 ± 0.33−52.86 ± 0.25
S380I−71.96 ± 0.41−40.69 ± 0.31
D421I−62.95 ± 0.39−45.20 ± 0.29
I389D−70.14 ± 0.31−47.45 ± 0.23
I389W−71.06 ± 0.34−45.55 ± 0.22
I389T−57.75 ± 0.39−42.92 ± 0.26
D454W−62.27 ± 0.37−44.12 ± 0.26
D421W−42.71 ± 0.37−41.08 ± 0.26
E40/ERK2−75.64 ± 0.27−49.50 ± 0.2
Table
Table 3. Comparison of the energetics in kcal/mol (±SEM) of mutants with E40/ERK2 obtained by the MM-PBSA/GBSA method. VDWAALS, van der Waals energy; EEL, MM electrostatic energy. The polar and non-polar contributions are EGB (or EPB) and ESURF (or ENPOLAR), respectively, for MM-GBSA (or MM-PBSA).
Table 3. Comparison of the energetics in kcal/mol (±SEM) of mutants with E40/ERK2 obtained by the MM-PBSA/GBSA method. VDWAALS, van der Waals energy; EEL, MM electrostatic energy. The polar and non-polar contributions are EGB (or EPB) and ESURF (or ENPOLAR), respectively, for MM-GBSA (or MM-PBSA).
MutantsA443DA443NS380LS380IN422AN422IN422T
EnergeticsVDWAALS−97.75 ± 0.23−88.88 ± 0.22−91.76 ± 0.25−97.73 ± 0.22−94.84 ± 0.28−86.68 ± 0.26−91.07 ± 0.24
EEL−626.46 ± 2.63−589.4 ± 2.97−673.80 ± 1.87−521.631 ± 2.08−441.77 ± 1.77−457.6 9 ± 1.88−407.50 ± 2.08
EGB672.13 ± 2.41634.15 ± 2.86712.20 ± 1.78582.05 ± 1.95503.72 ± 1.76512.63 ± 1.82448.52 ± 1.97
ESURF−7.78 ± 0.03−7.43 ± 3.01−3.38 ± 0.03−3.38 ± 0.04−5.97 ± 0.03−5.75 ± 0.04−2.81 ± 0.03
∆Ggas−724.22 ± 2.59−678.31 ± 3.01−765.56 ± 1.89−619.36 ± 2.07−536.62 ± 1.84−544.37 ± 1.92-498.58 ± 2.07
∆Gsolvation664.35 ± 2.41626.72 ± 2.86708.82 ± 1.77578.67 ± 1.94497.75 ± 1.75506.88 ±1.80445.71 ± 1.97
∆G TOTAL−59.8 6 ± 0.29−51.59 ± 0.26−56.74 ± 0.27−40.69 ± 0.31−38.87 ± 0.27−37.49 ± 0.27−52.86 ± 0.25
MutantsD421II389DI389WI389TD454WD421W
EnergeticsVDWAALS−102.03 ± 0.2−92.03 ± 0.22−93.21 ± 0.22−84.65 ± 0.22−104.74 ± 0.27−104.74 ± 0.27
EEL−357.28 ± 2.08−409.02 ± 2.03−501.96 ± 1.85−341.38 ± 2.02−307.17 ± 2.28−274.95 ± 1.61
EGB423.76 ± 1.81463.06 ± 1.95558.67 ± 1.79390.69 ± 1.95379.18 ± 2.29348.28 ± 1.54
ESURF−9.65 ± 0.01−9.04 ± 0.01−9.04 ± 0.01−8.77 ± 0.01−10.18 ± 0.01−10.06 ± 0.01
∆Ggas−459.31 ± 1.96−501.06 ± 2.01−595.18 ± 1.86−426.04 ± 1.97−411.92 ± 2.40−379.30 ± 1.62
∆Gsolvation414.11 ± 1.81453.61 ± 1.95549.62 ± 1.79381.91 ± 1.95368.99 ± 2.28338.21 ± 1.53
∆G TOTAL−45.20 ± 0.29−47.45 ± 0.23−45.55 ± 0.22−44.12 ± 0.26−42.92 ± 0.2−41.08 ± 0.26
Table
Table 4. Distance, occupancy, and angle of the hydrogen bonds formed from ERK2-A443D/D421W/S380L.
Table 4. Distance, occupancy, and angle of the hydrogen bonds formed from ERK2-A443D/D421W/S380L.
SystemAcceptorDonorOccupancy (%)DistanceAngle
A443D/ERK2ASP475-OD2TYR222-OH99.02.66166.77
ASP442-OD2HIE221-NE284.62.80154.14
ASP454-OD1ARG180-NH152.52.78159.34
VAL177-OARG455-NH150.72.81154.93
ASP224-OD2TYR444-OH43.92.70161.74
ASN488-OTYR176-OH42.72.76160.74
ASP409-OD2LYS220-NZ22.22.79152.82
TYR176-OARG455-NH16.22.86160.46
ASN190-OD1LYS479-NZ3.62.85154.21
D421W/ERK2ASP475-OD2TYR222-OH99.02.66166.77
ASP442-OD2HIE221-NE284.62.81154.13
VAL177-OARG455-NH150.72.81154.92
ASP454-OD2ARG180-NH146.72.80159.39
ASP224-OD1TYR444-OH43.92.71161.74
ASN488-OTYR176-OH42.72.75163.75
ASP454-OD1ARG180-NH138.32.80157.34
ASP409-OD2LYS220-NZ22.22.79152.82
TYR176-OARG455-NH16.22.85160.46
S380L/ERK2ASP475-OD2TYR222-OH98.92.65167.16
ASP442-OD2HIE221-NE284.42.80153.80
ASP421-OD1ARG180-NH272.62.80154.62
ASN488-OTYR176-OH68.72.75162.41
ASP454-OD1ARG180-NH155.12.77160.96
ASP224-OD1TYR444-OH52.82.70161.61
ASP421-OD1ARG180NH149.62.84149.43
TYR176-OARG455-NH145.32.85159.24
ASP409-OD2LYS220-NZ14.02.79154.31
Table
Table 5. Type of receptor-ligand interactions in mutant complexes.
Table 5. Type of receptor-ligand interactions in mutant complexes.
SystemAcceptorDonorInteraction
A443D/ERK2ASP475-OD2TYR222-OHCation-Π
ASP442-OD2HIE221-NE2Salt Bridge
ASP454-OD1ARG180-NH1Salt Bridge
VAL177-OARG455-NH1Salt Bridge
ASP224-OD2TYR444-OHCation-Π
ASN488-OTYR176-OHAmino-Π
ASP421-OD1ARG455-NHSalt Bridge
ASP409-OD2LYS220-NZSalt Bridge
TYR176-OARG455-NHCation-Π
ASN190-OD1LYS479-NZSalt Bridge
D421W/ERK2ASP475-OD2TYR222-OHCation-Π
ASP442-OD2HIE221-NE2Salt Bridge
VAL177-OARG455-NH1Salt Bridge
ASP454-OD2ARG180-NH1Salt Bridge
ASP224-OD1TYR444-OHCation-Π
ASN488-OTYR176-OHAmino-Π
ASP454-OD1ARG180-NH1Salt Bridge
ASP409-OD2LYS220-NZSalt Bridge
TYR176-OARG455-NHCation-Π
TRP421-OARG455-NECation-Π
S380L/ERK2ASP475-OD2TYR222-OHCation-Π
ASP442-OD2HIE221-NE2Salt Bridge
ASP421-OD1ARG180-NH2Salt Bridge
ASN488-OTYR176-OHAmino-Π
ASP454-OD1ARG180-NH1Salt Bridge
ASP224-OD1TYR444-OHCation-Π
ASP421-OD1ARG180NH1Salt Bridge
TYR176-OARG455-NH1Cation-Π
ASP409-OD2LYS220-NZ2Salt Bridge
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Gautam, V.; Nimmanpipug, P.; Zain, S.M.; Rahman, N.A.; Lee, V.S. Molecular Dynamics Simulations in Designing DARPins as Phosphorylation-Specific Protein Binders of ERK2. Molecules 2021, 26, 4540. https://doi.org/10.3390/molecules26154540

AMA Style

Gautam V, Nimmanpipug P, Zain SM, Rahman NA, Lee VS. Molecular Dynamics Simulations in Designing DARPins as Phosphorylation-Specific Protein Binders of ERK2. Molecules. 2021; 26(15):4540. https://doi.org/10.3390/molecules26154540

Chicago/Turabian Style

Gautam, Vertika, Piyarat Nimmanpipug, Sharifuddin Md Zain, Noorsaadah Abd Rahman, and Vannajan Sanghiran Lee. 2021. "Molecular Dynamics Simulations in Designing DARPins as Phosphorylation-Specific Protein Binders of ERK2" Molecules 26, no. 15: 4540. https://doi.org/10.3390/molecules26154540

Article Metrics

Back to TopTop