Abstract
Background
Extracellular signal-regulated kinase 1 (ERK1) belongs to mitogen-activated protein kinases, which are essential for memory formation, cognitive function, and synaptic plasticity. During Alzheimer's disease (AD), ERK1 phosphorylates tau at 15 phosphorylation sites, leading to the formation of neurofibrillary tangles. The overactivation of ERK1 in microglia promotes the release of pro-inflammatory cytokines, which results in neuroinflammation. Additionally, elevated oxidative stress during AD stimulates the ERK1 pathway, leading to neuronal loss.
Objective
Because ERK1 signaling plays a significant role in tau phosphorylation, targeting ERK1 may be therapeutically beneficial by either preventing excessive activation of the signaling pathway or altering its pathway to enhance neuroprotective effects during AD.
Methods
This study employed structure-based virtual screening of phytoconstituents from the IMPPAT library. Subsequently, in-depth docking and molecular dynamics (MD) simulation studies were implemented to identify potential ERK1 inhibitors with desirable pharmacological properties.
Results
Silandrin and Hydroxytuberosone were found to be potential ERK1 inhibitors with higher affinity and specificity than the control molecule Tizaterkib. These compounds specifically bind to the ERK1 substrate binding pocket and interact with crucial residues. Finally, the elucidated compounds with ERK1 were evaluated using an all-atom molecular MD simulation to analyze structural dynamics, structural compactness, hydrogen bond dynamics, principal component analysis, and free energy landscape.
Conclusions
The study suggested that Silandrin and Hydroxytuberosone can further be exploited as potential lead molecules for therapeutic development against ERK1-mediated AD.
Keywords
Introduction
Alzheimer's disease (AD) is a neurodegenerative disorder characterized by the progressive decay of neuronal cells, leading to a gradual worsening of dementia symptoms. Dementia involves the loss of cognitive abilities such as memory, language, and reasoning, severely impacting the daily lives of affected individuals.1,2 AD and various forms of dementia affect approximately 55 million people globally, with AD accounting for 60–80% of all dementia cases (http://www.alz.org). In many instances, dementia is linked to multiple underlying changes in the brain. 3 The number of individuals affected by AD is expected to double every five years, reaching an estimated 152 million by 2050. 4
The research community has proposed several hypotheses to explain the potential causes of AD, including the cholinergic hypothesis, amyloid hypothesis, and tau hypothesis, among others. This study focuses on the tau hypothesis, which posits that the abnormal function of the microtubule-associated tau protein, specifically tau hyperphosphorylation, is responsible for the onset of AD. 5 Tau hyperphosphorylation leads to the aggregation of tau protein into neurofibrillary tangles (NFTs) by forming paired helical filaments of tau protein. NFTs are a defining characteristic of AD.6,7 In normal neurophysiological conditions, the microtubule maintains the cytoskeleton of neurons in the brain, and the tau protein regulates the function of the microtubule. For ideal protein activity, tau protein has 30 phosphorylation sites at which phosphorylation occurs. 8 In affected individuals, due to hyperphosphorylation of tau-protein, disruption in the regulation of microtubules takes place, and it results in the decay of neuronal cells.9,10
Extracellular signal-related kinase 1 (ERK1), also known as mitogen-activated protein kinase 3 (MAPK3), is a promising drug target because its hyperactivity causes tau-hyperphosphorylation. 11 It belongs to the superfamily mitogen-activated protein kinase (MAPK). It is a proline-directed kinase that phosphorylates serine/threonine residues, followed by a proline residue (Ser/Thr-Pro) motif of tau-protein. This protein kinase is encoded in gene ERK1 on chromosome 16. It is localized in the cytoplasm, nucleus, membrane caveola, cell junction, and focal adhesion. 12 ERK1 consists of 379 amino acid residues. The protein kinase domain spans residues Y42 to L330. The adenosine triphosphate (ATP) binding site includes residues I48 to V56 and K71. The active site, which attracts a proton, is located at residue D166. The TXY motif, essential for the activation of MAP kinases, is found between residues T202 and Y204. The phosphorylation of the threonine and tyrosine residues within the TXY motif is crucial for MAP kinase activation. The cofactor for ERK1 is Mg2+ ion.13,14
There is no definitive cure for the treatment of dementia at present. 15 Currently, there are symptomatic treatments available for AD, classified into two categories of drugs: cholinesterase inhibitors and N-methyl-D-aspartate (NMDA) receptor antagonists. These medications help manage the symptoms associated with AD but do not cure the disease. 16 In this study, disease-modifying therapy (DMT) was employed to combat the progression of tau-associated AD. Unlike symptomatic treatments, DMT induces significant changes in the clinical progression of AD by targeting the underlying pathophysiological pathways that lead to neuronal cell decay. 17 Targeting the hyperactivity of ERK1 is anticipated to be a promising DMT for AD. Activation of ERK1 contributes to the phosphorylation of tau, leading to its hyperphosphorylation, ultimately leading to the development of NFTs.18,19 Inhibiting ERK1 can resist the hyperphosphorylation of tau, which in turn will reduce the chances of AD progression.
This approach aims to identify novel therapeutic leads from the Indian Medicinal Plants, Phytochemistry, and Therapeutics (IMPPAT 2.0) library that inhibit ERK1, which is involved in the tau-associated pathology of AD, thereby modifying the disease and potentially providing a definitive cure for tau-associated AD progression—a significant challenge in the past. Natural products offer various advantages in targeted drug discovery. They show structural diversity, desired biological activity, and less toxicity than synthetic compounds. Natural products are being evolutionary optimized to interact with biological macromolecules with high efficacy and selectivity as therapeutic agents. The small molecules derived from plant sources, as well as other natural sources, have also shown historical success in targeted drug discovery. 20
Medicinal plants have been used for centuries in various traditional medical systems to treat neurological diseases.21,22 The majority of identified natural products are either nontoxic or therapeutic, but some have been shown to have harmful effects. As a result, the source and degree of toxicity of natural products must be carefully considered. An expanded and enhanced phytochemical atlas of Indian medicinal plants, IMPPAT 2.0, offers a non-redundant, FAIR24-compliant, in silico stereo-enabled library of 17,967 phytochemicals. Indian medicinal plants have long been used in traditional Indian medical systems like Ayurveda and Siddha to treat a range of diseases.
This research employs a structure-based virtual screening approach to identify novel inhibitors with therapeutic potential against ERK1, as illustrated in Figure 1. Virtual screening, coupled with molecular dynamics (MD) simulation, has significantly streamlined the process of identifying potential inhibitors with desirable pharmacological properties from extensive libraries. MD simulation is particularly valuable as it offers critical insights into the dynamic behavior of the kinase-inhibitor complex at the atomic level and the energetics of complex interactions, ultimately validating potential inhibitors against the targeted kinase.

Workflow utilized in the identification of potential inhibitors against ERK1.
Methods
Computer environment and the online resources
This study was conducted on Asus TUF gaming F15, a local workstation running Windows 11 with a 2.5 GHz processor, 8 GB of RAM, 4GB of dedicated graphics processing unit (GPU), and a 512 GB SSD. A dependable, steady power supply was used while using high-speed internet. Tools like PyMOL 2.5 23 and InstaDock 1.0 24 were used to visualize the protein structure, preparation, and molecular docking of protein-ligand complexes, respectively. A workstation running on the Linux operating system with 24 GB of dedicated GPU and 32 GB of RAM was used for MD simulation. For the retrieval, assessment, and analysis of hit molecules, this study used various online resources, such as UniProt for analysis of research done on AD, RCSB-PDB for collecting the 3D structure of ERK1, 25 IMPPAT 2.0 database for collection of phytochemicals, 26 SwissADME for prediction of medicinal chemistry, drug-like properties, water solubility etc. of hit molecules, 27 pkCSM for prediction of ADMET properties of hit molecules, 28 Carcinopred-EL for prediction of carcinogenicity of hit molecules, 29 way-2-drug for prediction of activity spectra for hit molecules. 30 Discovery Studio 2024 software was employed for the interaction analysis of protein-hit molecules complex. 31
Receptor preparation
Using UniProt and the Protein Data Bank (PDB), the appropriate ERK1 structure was selected based on high resolution, the presence of the kinase domain, associated ligands, and structural purity. The human ERK1 X-ray crystal structure (PDB ID: 4QTB, Resolution: 1.40 Å) was downloaded from the RCSB PDB. The structure was further refined using PyMOL 2.5, where water molecules and heteroatoms, including co-crystallized ligands, were removed. Since the 3D structure of ERK1 was complete and free of gaps, homology modeling was unnecessary. A 3D configuration grid encompassing the entire protein was generated for binding hit molecules using AutoDock Tools 1.5.7. The prepared receptor file was then saved in the Protein Data Bank with partial charge and atom type (PDBQT) format and subjected to molecular docking using InstaDock.
Compound library used for virtual screening
A local library was assembled for plant-based hit molecules, comprising 3D structures of 11,708 molecules in ‘.mol’ file format sourced from the IMPPAT 2.0 database. To ensure the selection of suitable compounds, a control compound (AZD0364) was obtained from PubChem using an advanced search algorithm based on the Lipinski Rule of Five, a fundamental protocol in small molecule therapeutics. The Lipinski Rule delineates specific parameters and threshold values that characterize drug-like compounds. According to this rule, compounds should possess a molecular weight below 500 Da, with no more than 5 hydrogen bond donors, 10 hydrogen bond acceptors, and a LogP value below 5. These criteria help identify molecules with favorable pharmacokinetic properties for further investigation. Tizaterkib, which is a known inhibitor of ERK1, is used as a reference molecule in this study to validate the findings. The selection of Tizaterkib as the reference compound for ERK1 inhibition was based on several factors that make it suitable for comparison with the identified natural product inhibitors, Silandrin and Hydroxytuberosone. Tizaterkib is a clinically relevant ERK1 inhibitor that has shown potential in inhibiting the MAPK pathway in cancer research. 32
Tools and directories used in molecular docking
We employed InstaDock 1.0 for molecular docking-based virtual screening of hit molecules within the local library. InstaDock serves as a standalone software platform offering a comprehensive suite for molecular docking. Utilizing OpenBabel 3.0 and RasMol packages in the backend facilitates file conversion and visualization, respectively. A designated working directory is essential to initiate molecular docking for virtual screening. This directory stores the InstaDock, the receptor file, and the ligands. Upon program execution, InstaDock first generates a conf file to establish a 3D space for ligand docking. Subsequently, it converts receptor and ligand files into the requisite receptor and ligands in ‘.pdbqt’ format for molecular docking. Upon completion of the docking process, InstaDock generates a “Result” directory within the working directory. This directory encompasses several files for each docked ligand, including ligand_log, ligand_out, and affinity_result. Additionally, it creates a subdirectory named “top_hits,” housing further analysis such as affinity_top_hits, ligand_log, ligand_out, and top_split_files for each ligand.
Molecular docking-based virtual screening
The 3D structure of the target ERK1 (PDB-ID: 4QTB) was retrieved from the Protein Data Bank, followed by virtual screening using the receptor-ligand molecular docking approach. Molecular docking enabled the observation of ligand orientation within the protein's active site, shedding light on the mechanism underlying structure-based drug design and substrate inhibitor selectivity. A library of phytochemicals with drug-like properties sourced from the IMPPAT 2.0 database in 3D (.mol) format was utilized alongside the 3D structure of a control compound obtained from PubChem. The ligand preparation was performed using InstaDock. The ligand preparation steps involve the removal of the non-polar hydrogen atoms and the addition of Gasteiger charges. Then, the root of the ligand's torsion tree was detected, and the number of torsions was chosen. Finally, the prepared ligand file is saved in PDBQT format. Molecular docking facilitated the examination of bond conformations and ligand binding affinities with ERK1. To establish optimal ligand-receptor complexes, ligands were docked using a predefined grid box with specific 3D configuration parameters (x center = 13.233, y center = 9.614, z center = 16.671) and size parameters (x = 69, y = 70, z = 63). Various ligand orientations were explored within the kinase domain of the protein, emphasizing the ATP binding site (residues I48-V56, K71) and the active site (residue D166) of the receptor. The docking algorithm employed generated a random seed for each ligand screening, producing a scoring function to assess ligand binding affinity. Following virtual screening completion, docking results were analyzed from the generated .out and .log files, identifying the most suitable docked conformations based on high binding affinity for further investigation. For visualization and interaction analysis of the receptor-ligand complex, PyMOL 2.5 and Discovery Studio Visualizer 2024 were utilized, enabling a comprehensive exploration of molecular interactions and structural insights.
Analysis of medicinal chemistry, physicochemical properties
After the molecular docking, molecules displaying significant affinity scores were identified and subsequently subjected to comprehensive analysis using the SwissADME web-based tool. This analysis encompassed an examination of medicinal chemistry, physicochemical attributes, water solubility, and drug-like properties, crucial for understanding the fate of hit molecules within organisms. Medicinal chemistry analysis included the prediction of PAINS (Pan Assay Interference Compounds), ensuring the exclusion of compounds prone to interference in assays. Physicochemical properties were assessed, with a particular focus on predicting polar surface area utilizing the technique of topological polar surface area (TPSA). Evaluation of drug-likeness was conducted based on Lipinski's rule, a fundamental criterion in assessing compound viability for drug development. Additionally, an analysis of water solubility was performed, considering parameters such as the ESol class, Ali class, and Silicos-IT class, to gauge the compound's potential solubility in biological systems.
ADMET analysis
After the assessment of physicochemical, drug-likeness, and water-solubility properties, compounds exhibiting favorable scores for desired attributes were further scrutinized through absorption, distribution, metabolism, excretion, and toxicity (ADMET) analysis. ADMET was conducted using the web-based tool pkCSM. 33 Compounds demonstrating significant scores for desired properties related to ADMET while exhibiting minimal toxicity concerns such as hepatotoxicity, Salmonella typhimurium reverse mutation assay (AMES) toxicity, and skin sensitivity were prioritized for subsequent analysis.
Analysis of carcinogenicity and biological activity
After ADMET analysis, hit molecules with a significant score for ADMET properties were analyzed for carcinogenic activities using the web-based tool CarcinoPred-EL. The study was based on the threshold value for chemical carcinogenicity, i.e., (average probability = 0.5). If the average probability of a compound is greater than 0.5, then it is carcinogenic and vice-versa. Non-carcinogen compounds were analyzed for their biological activity using the web-based tool Way-to-drug. Compounds showing therapeutic activities were selected from all biological activity predicted in PASS analysis. The compounds were selected based on Pa (probability of being active) and Pi (probability of being inactive).
Interaction analysis
All hit molecules having non-carcinogenic activity were subjected to interaction analysis. An interaction analysis was performed using Discovery Studio software. All 9 conformers of each compound, along with the control, were imported into the discovery studio along with the .pdb file of the receptor protein. Receptor-ligand complexes were analyzed, and complexes with more hydrogen bonds in the active or binding site were selected. 2D diagrams were generated to investigate the interacting bonds formed when the ligand interacts with the receptor. Complex interactions with hydrogen bonds were considered more stable than other non-covalent bonds.
Molecular dynamic simulations
MD simulations are a crucial tool in computer-aided drug development for studying the movements and interactions at the atomic level in receptor-ligand complexes.34–37 MD simulations validated the stability of the docked complex of ERK1-IMPHY011108, ERK1-AZD0364, and ERK1-IMPHY014037 generated after molecular docking. GROMACS 2022.2 package was used for MD simulations. PRODRG server was used to create the topology of ERK1, ERK1-AZD0364, ERK1-IMPHY011108, and ERK1-IMPHY014037 complexes. 38 The solvation of each system was achieved by placing them in a cubic box, where the distance from the center to the edges was set to 10 Å. The simple point charge (SPC216) water model was used for the solvation process, and 26,142 water molecules were added to the system for simulation. In addition, the simulated system was made neutral by introducing 3 Na+ counter-ions.
An energy minimization procedure was carried out to alleviate any potential clashes between atoms. The solvated system underwent 1500 steps of energy minimization using the steepest descent method, followed by the conjugate gradient method. For equilibration, a two-step process was employed by utilizing periodic boundary setting. In two steps, the system was initially equilibrated at a constant volume for 1000 ps. During the first time, the temperature was gradually increased from (0–300 K), and the second equilibration was performed by maintaining a pressure of 1 atm. At last, an MD simulation was performed for 100 ns.
Results
Molecular docking-based virtual screening
Using a molecular docking approach, we conducted a virtual screening of phytochemicals alongside a control compound to identify competitive inhibitors of ERK1. The information about the binding site and active site of ERK1 has been extracted from the UniProt database. The ATP-binding site of the ERK1 protein is located between the residues Gln46-Ser58 and Lys71. Additionally, the active binding site of ERK1 is situated at Asp166. 39 Our screening revealed numerous phytochemicals with robust binding affinity scores toward the ERK1 binding pocket, suggesting their potential as ERK1 inhibitors. Among the 11,708 compounds screened, we identified 80 top hits. These hits exhibit binding affinity values ranging from −12.1 to −9.3 with ERK1, indicating significant potential for interaction with this target (Supplemental Table 1).
These compounds displayed appreciable binding affinity, warranting further investigation to ascertain their inhibitory potential against ERK1. These 80 hits were filtered for their physicochemical properties, where 18 compounds showed appreciable drug likeliness with zero PAINS, as discussed in the ensuing section. In Table 1, we present the binding affinities and sources of the top 18 hits that offer insight into the promising candidates emerging from docking screening.
Binding affinities of top 18 hits and control compound binding with ERK1.
Analysis of physicochemical properties
The hits generated from the molecular docking screening were evaluated for their physicochemical properties using SwissADME. Here, a total of 62 compounds were eliminated based on their physicochemical properties (TPSA-topological polar surface area), water solubility (log S-ESOL, log S-Ali, log S-SILICOS-IT), drug-likeness (Lipinski violation), medicinal chemistry (pan-assay interference compounds structure - PAINS). The threshold value for TPSA was set to (65–145), water solubility was set to above poor solubility, Lipinski violation was set to 0, and PAINS was set to 0. As compounds with PAINS structure can lead to undesired regulation and outcomes, compounds containing PAINS were eliminated. Predicted SwissADME parameters of 18 screened compounds are shown in Table 2.
Properties of hit molecules predicted by SwissADME.
TPSA: topological polar surface area.
ADMET analysis
ADMET analysis predicts the pharmacokinetic properties of hit molecules.40,41 Phytochemicals selected from the initial physicochemical screening were subjected to ADMET analysis using the pkCSM server. We have analyzed the human intestinal absorption and blood-brain barrier permeation for the compounds. For metabolism, CYP2D6 inhibition was checked, and for excretion, OCT2 substrate was accounted for. Further, we selected those compounds that showed negative results for hepatotoxicity and AMES toxicity. Out of the evaluated compounds, nine successfully passed the ADMET screening criteria. All assessed parameters fell within acceptable ranges for potential drug candidates, as detailed in Table 3. This screening ensures that these compounds exhibit favorable pharmacokinetic properties, which are crucial for further drug development processes.
ADMET properties of hit molecules predicted by pkCSM.
Analysis of carcinogenicity and biological activity
None of the compounds screened from the ADMET analysis further underwent CarcinoPred-EL screening to predict any carcinogenic patterns. The study found all compounds free of any carcinogenic patterns. The threshold for chemical carcinogenicity is set at an average probability of 0.5. The average probability across all 9 compounds exceeded this threshold. Detailed average probabilities for the screened compounds are presented in Table 4. These non-carcinogen compounds were subjected to the prediction of activity spectra for substances (PASS) analysis by utilizing a web-based tool, way-to-drug, to predict biological activity by importing the SMILES of each compound. PASS analysis gave a set of biological activities based on Pa and Pi values. 42 As ERK1 is responsible for the hyperphosphorylation of tau protein, the proliferation of neuronal cells, and causing dementia, therapeutic activities were chosen antagonists to activity shown by receptors in the affected individual. Here, three compounds (IMPHY011108, IMPHY014037, and IMPHY014814) were selected with considerable therapeutic activities (Table 5). The PASS results showed that these compounds can be exploited in the therapeutic development against ERK1-associated cancer.
The average probability of carcinogenicity of hit molecules predicted by Carcinopred-EL.
The biological activity of hit molecules is computed by PASS analysis utilizing Way2Drug.
Interaction analysis
The screened phytochemicals and the control compound with considerable therapeutic activities were selected for interaction analysis. All possible docked conformations for each compound and their interaction with ERK1 were analyzed using PyMOL and Discovery Studio Visualizer software. Here, two compounds, IMPHY011108 (Silandrin) and IMPHY014037 (Hydroxytuberosone), were selected based on their specific interaction with the critical residues of active and binding sites of ERK1. The active site for ERK1 is located at residue D166, while ATP binding sites are found between residues I48-V56 and at residue K71. The binding of hit molecules to these ATP binding sites on ERK1 potentially inhibits ATP availability, thereby inhibiting ERK1 activity. The selection of the best binding pose of the ligands was based on a combination of molecular binding affinity and the specific interaction patterns with key residues of ERK1. The selected binding pose should establish hydrogen bonds or other strong non-covalent interactions with these residues, which are essential for the effective inhibition of ERK1 activity. Using the PyMOL tool, we visualized the complexes of Silandrin, Hydroxytuberosone, and AZD0364 (Tizaterkib) with ERK1 and generated a vacuum electrostatic field for these complexes. A detailed view of the 3D interactions and vacuum electrostatic fields of all complexes with ERK1 is presented in Figure 2.

Protein-ligand interactions. (A) 3D-Structural representation of ERK1-kinase domain in complex with Silandrin (red), Hydroxytuberosone(violet), and Tizaterkib (green). (B) Detailed representation of protein-ligand complex showing interacting residues. (C) Electrostatic potential representation of ERK1 with its complexes.
The two-dimensional interactions of Silandrin, Hydroxytuberosone, and Tizaterkib with the interacting residues of ERK1 were plotted and are shown in Figure 3. The plots showed different types of interactions between the protein and ligands, such as Van der Waals forces, conventional hydrogen bonds, carbon-hydrogen bonds, pi-anion interactions, pi-sigma interactions, pi-pi T-shaped interactions, alkyl interactions, halogen (fluorine) interactions, and pi-alkyl interactions. Silandrin established conventional hydrogen bonds with the ATP binding site residues Y53 and K71 (Figure 3A). Hydroxytuberosone formed conventional hydrogen bonds with residues Y53, N171, and S170 (Figure 3B). Tizaterkib established a traditional bond of hydrogen with residue D184 (Figure 3C). Silandrin forms a pi-sigma bond with Leu173 and a pi-anion bond with Glu88. Hydroxytuberosone forms a pi-sigma bond with Val56, which is located in the ATP binding region.

Structural representation in 2D depicting ERK1 interaction with (A) Silandrin, (B) Hydroxytuberosone, and (C) tizaterkib. The complexes show different interactions with the residues of the binding pocket of ERK1, where Tyr53 and Lys71 are ATP binding sites.
Molecular dynamic simulations
MD simulation is a useful procedure in computer-aided drug discovery for analyzing physical interaction in a biophysical environment based on structural dynamics, structural compactness, hydrogen bond dynamics, principal component analysis (PCA), and analysis of free energy landscape (FEL) for the complexes of phytochemicals and receptors. 43 After interaction analysis, the elucidated compounds Silandrin and Hydroxytuberosone showed stable interactions on binding pockets of ERK1. These compounds, along with the control molecule Tizaterkib, were further validated utilizing the MD simulation approach in an explicit solvent for its physical interaction as a function of time for 100 ns.
Structural dynamics
The protein structural dynamics were investigated before and after the phytochemicals were bound to ERK1. The structural dynamics of protein were carried out by analyzing two parameters, root-mean-square deviation (RMSD) and root-mean-square fluctuation (RMSF). RMSD is a prominent statistic for evaluating protein structural deviation. It is one of the most important parameters for analyzing the deviation of a protein's backbone from its native conformation at the level of atomic coordinates. 44 ERK1, complexes of ERK1-Silandrin, ERK1-Hydroxytuberosone, and ERK1-Tizaterkib were analyzed, and it was found that the mean value of probability distribution function (PDF) for RMSD was 0.33 nm, 0.38 nm, 0.25 nm, and 0.27 nm, respectively. The binding of Silandrin to ERK1 showed the least deviation from ERK1, which indicated that the docked complex was robust and at equilibrium throughout the simulation, as evidenced by the RMSD graph in Figure 4A.

Time-evolution plots demonstrating structural dynamics of ERK1 and its complex with Silandrin, Hydroxytuberosone, and Tizaterkib. (A) Root-mean-square-deviation (RMSD) of ERK1 and its complexes. (B) The probability distribution function (PDF) of RMSD. (C) Residual fluctuation (RMSF) of ERK1 and its docked complexes. (D) The PDF of RMSF values.
RMSF is another useful parameter for assessing the residual vibrations in a protein structure during an MD simulation, as it helps get insights into the effects of ligand binding on the residual fluctuations of the protein. It sheds light on the residual dynamics of ERK1 before and after ligand binding by plotting the RMSF for each residue. ERK1, ERK1-Silandrin, ERK1-Hydroxytuberosone, and ERK1-Tizaterkib were analyzed for their RMSF, where it was found that the mean value of PDF for RMSF was 0.14 nm, 0.14 nm, 0.12 nm, and 0.15 nm respectively. The complex of ERK1 with Silandrin and Hydroxytuberosone showed identical residual fluctuation patterns as in free ERK1. In contrast, Tizaterkib, upon binding with ERK1, shows a considerable difference in residual fluctuation in comparison with free ERK1, which signifies that Silandrin and Hydroxytuberosone bind ERK1 more robustly (Figure 4B).
Structural compactness
The structural compactness of ERK1 and its complexes were studied in time-evolution settings before and after the binding of compounds. Structural compactness was studied based on two parameters: Radius of gyration (Rg) and solvent-accessible surface area (SASA). Rg is a valuable measure for analyzing protein structural compactness based on the orientation of atoms along the central axis of the protein. 45 It is directly related to the protein structure's folding and compaction. Rg is the root mean square (RMS) distance between the set of atoms and their collective center of mass. After analyzing the PDF of Rg for ERK1 and its complexes with Silandrin, Hydroxytuberosone, and Tizaterkib, it was found that the mean value of Rg was 2.20 nm, 2.22 nm, 2.19 nm, 2.15 nm, respectively. ERK1- Hydroxytuberosone shows the least difference from ERK1, which indicates that the binding of Hydroxytuberosone with ERK1 shows stable structural compactness of its complex. In terms of Rg, the binding of Silandrin with ERK1 is less stable than Hydroxytuberosone, and Tizaterkib has shown considerable difference (Figure 5A).

Time-evolution plots demonstrating the structural compactness of ERK1 and its complex with Silandrin, Hydroxytuberosone, and Tizaterkib. (A) The radius of gyration (Rg). (B) The probability distribution function (PDF) of Rg values. (C) Solvent accessible surface area (SASA). (D) The PDF of SASA values.
SASA is a protein's surface area accessible to its nearby solvent. SASA is a protein stability and folding behavior assay. When ERK1 and its complexes with Silandrin and Hydroxytuberosone, Tizaterkib were analyzed using SASA, the mean values of PDF found were 180.4 nm2, 181.6 nm2, 180.7 nm2, and 178.1 nm2, respectively. It was found that both complexes with Silandrin and Hydroxytuberosone were almost in equilibrium with the free ERK1 as they have a small difference in their SASA values; in contrast, Tizaterkib has shown a considerable difference, which is illustrated in Figure 5B. Out of both complexes, ERK1-Hydroxytuberosone is more proximal to free ERK1, and this justifies that Hydroxytuberosone can form stable interaction based on the SASA. Based on both parameters, i.e., Rg and SASA, Hydroxytuberosone binds with ERK1, and its complex is stable regarding overall structural compactness.
Hydrogen bond dynamics
Hydrogen bonds between the protein's backbone and amide nitrogen provide stability to the secondary structure of the protein. The structural strength of the protein's backbone is linked to hydrogen bonding dynamics. Using the MD simulation approach, the average hydrogen bond between ERK1 and compounds was estimated for whole trajectories. Two types of hydrogen bonds are present in complex, i.e., intramolecular hydrogen bonds and intermolecular hydrogen bonds. Intramolecular hydrogen bonds are essential in regulating the integrity of a protein's structure. 46 The intramolecular hydrogen bonds were analyzed over a 100 ns simulation period for ERK1 and its complexes with Silandrin, Hydroxytuberosone, and Tizaterkib (Figure 6). The average number of intramolecular hydrogen bonds observed was 268 for ERK1, 266 for the ERK1-Silandrin complex, 271 for the ERK1-Hydroxytuberosone complex, and 276 for the ERK1-Tizaterkib complex. Among these compounds, Silandrin exhibited the smallest deviation from the ERK1 baseline, indicating the most stable dynamics of intramolecular hydrogen bonds (Table 6).

Plot showing theintramolecular hydrogen bonds for ERK1 and its complex with Silandrin, Hydroxytuberosone, and Tizaterkib. (A) Time evolution plot of the number of intramolecular hydrogen bonds. (B) Probability distribution function plot of intramolecular hydrogen bonds.
Mean values of MD parameters computed for 100 ns simulation. # represents the intramolecular hydrogen bond number.
The dynamics of intermolecular hydrogen bonding play a crucial role in stabilizing ERK1 and its complexes with Silandrin, Hydroxytuberosone, and Tizaterkib. Throughout the simulation, these hydrogen bonds act as stabilizing forces, maintaining the integrity of the docked complexes and preventing significant deviations from the initial docked positions. Upon detailed analysis, it was found that the probability density function (PDF) values for intermolecular hydrogen bonds were higher in the ERK1 complexes with Hydroxytuberosone and Silandrin, indicating stronger and more frequent interactions. This suggests that these compounds form more stable complexes with ERK1. In contrast, the ERK1-Tizaterkib complex exhibited the lowest PDF values for intermolecular hydrogen bonds, suggesting fewer and less stable interactions, as shown in Figure 7. These findings underscore the importance of intermolecular hydrogen bonding in drug design, highlighting how different compounds can vary significantly in their ability to stabilize target proteins. This information is crucial for optimizing lead compounds in drug development, as it provides insights into the molecular interactions that contribute to complex stability.

Time-evolution plots demonstrating analytics of intermolecular hydrogen bonds for ERK1 with (A) Silandrin, (B) Hydroxytuberosone, and (C) Tizaterkib. The lower panel shows the probability distribution function of intermolecular hydrogen bonds.
Principal component analysis
PCA is applied to analyze the dynamics of the collective motion for ERK1 and its complexes with Silandrin, Hydroxytuberosone, and Tizaterkib.47,48 Using trajectories of simulated ERK1 and its complexes, analysis of conformation formed on distinct eigenvectors (EV1 and EV2) for ERK1, ERK1-Silandrin, ERK1-Hydroxytuberosone, and ERK1-Tizaterkib was performed (Figure 8). The generated plot depicts that Tizaterkib has a similar conformation to hit molecules with the free conformation of ERK1 (Figures 8A and D). At the same time, ERK1-Silandrin and ERK1-Hydroxytuberosone show slightly different conformation space from ERK1 (Figure 8B and C). It can be anticipated that Hydroxytuberosone is more like ERK1 in conformational space in contrast with Silandrin; hence, the complex of ERK1 with Hydroxytuberosone is more stable.

Two-dimensional conformational projections of (A) ERK1 and its complex with (B) Silandrin, (C) Hydroxytuberosone, and (D) Tizaterkib computed from principal component analysis.
Analysis of free energy landscapes
This analysis provided insight into FELs characterizing the folding dynamics of ERK1 in a free state and its docked complexes with Silandrin, Hydroxytuberosone, and Tizaterkib. FELs provide insights into protein folding mechanisms by mapping the potential energy of different protein conformations, revealing energy barriers, pathways, and stable states.49,50 The FELs were represented as contoured maps, demonstrating distinct color gradients, with deeper blue regions representing energetically favorable conformations closely resembling the free state of ERK1 (Figure 9). In its free state, ERK1 exhibited a single dominant global minimum within a broad basin, which indicates its stability (Figure 9A).

3D-plot demonstrating Gibbs free energy landscape of (A) ERK1 and its complex with (B) Silandrin, (C) Hydroxytuberosone, (D) Tizaterkib.
However, upon binding to Silandrin, Hydroxytuberosone, and Tizaterkib, notable fluctuations were observed in both the size and positioning of the stable global minima (Figure 9B–D). The darker blue regions denoted regions of low energy and structural stability, suggesting that ERK1 remained predominantly in its native conformational state. Despite the introduction of Silandrin and Hydroxytuberosone, which induced the emergence of multiple basins and diverse conformations, the protein ultimately converged to its global minimum without undergoing significant unfolding. This analysis underscores ERK1's remarkable resilience in maintaining its native state even in the presence of interacting compounds. Hence, it can be anticipated that Silandrin and Hydroxytuberosone possess the potential to bind with ERK1 without inducing disruptive structural changes (Table 7). The conformational stability analysis parameters suggest a more stable complex formation of ERK1 with Hydroxytuberosone in terms of intramolecular hydrogen bonds. All these structural parameters indicate a more stable interaction in ERK1-Hydroxytuberosone. Overall, all the insights validated a promising foundation for Silandrin and Hydroxytuberosone as potential drug leads targeting ERK1.
Chemical and molecular properties of the finally selected phytochemicals against ERK1.
Discussion
AD presents a multifaceted challenge with various risk factors contributing to its onset and progression. DMTs targeting tau-hyperphosphorylation have been explored to intervene in AD progression. Among potential targets, ERK1 has emerged as pivotal due to its role in facilitating tau-hyperphosphorylation. Recent research on AD has demonstrated the importance of focusing therapeutic intervention on the ERK1 signaling pathway and tau protein. A lot of investigations have looked into ERK1 inhibition as a possible measure to manage tau pathology. Preclinical research has indicated that using particular inhibitors to target the ERK1 pathway may reduce tau phosphorylation and tau aggregation. 51 Also, drug candidates that modulate the ERK1 pathway have demonstrated neuroprotective effects in AD models, which showed improvement of cognitive function by reducing tau accumulation and promoting neuronal survival. 52
Taking the advantage of structural features of ERK1, structure-based virtual screening was conducted using phytochemicals from the IMPPAT 2.0 database. This approach efficiently narrowed down the phytochemical library, prioritizing molecules with potential ERK1 inhibitory activity. A library of 11,708 phytochemicals from the IMPPAT 2.0 database and a control compound, Tizaterkib (AZD0364), were retrieved for structure-based virtual screening against ERK1. These compounds were screened for various drug-like properties and binding potential towards ERK1. The results from each screening step are discussed in the ensuing sections.
There were a few inhibitors of ERK1, including Sulindac, 53 Purvalanol, Ulixertinib, 54 and Tizaterkib. We chose Tizaterkib as a reference molecule for our study. The reference molecule Tizaterkib showed a docking score of −9.5 kcal/mol. The most favorable conformation for each compound was identified by comparing their docked positions to that of Tizaterkib. The results showed that Silandrin and Hydroxytuberosone interact at the ATP binding sites of ERK1, specifically at residues Y53 and K71, with docking scores of −11.2 and −10.7 kcal/mol, respectively. These two compounds also interact with several other important residues of the ERK1 binding pocket. A detailed interaction analysis was carried out to understand the noncovalent interactions of Silandrin, Hydroxytuberosone, and the reference inhibitor Tizaterkib with ERK1.
In contrast, the control compound, Tizaterkib, exhibited notable deviations from native ERK1 dynamics. These findings indicate that Silandrin and Hydroxytuberosone could act as potent inhibitors of ERK1, potentially mitigating tau-associated AD progression. Silandrin, a flavonolignan extracted from Silybum marianum, and Hydroxytuberosone, derived from Pueraria tuberosa, exhibit therapeutic potential with demonstrated stability and binding affinity to ERK1.
The MD simulation study further established a correlation between the two ligands and the reference molecule Tizaterkib, showing comparable structural deviation, compactness, and FELs. The flavonolignan Silandrin has been shown in earlier research to have therapeutic potential for the treatment of neurodegenerative diseases. 55 Ancient medicine also used hydroxytuberosone as a bioactive ingredient to treat AD and other geriatric conditions. 56 Overall, the study suggests that Silandrin and Hydroxytuberosone hold significant potential as lead molecules for therapeutic development against ERK1-associated neurodegenerative diseases, including AD.
Conclusions
We have identified Silandrin and Hydroxytuberosone as promising candidates, both showing significant binding affinity to ERK1. Subsequent MD simulations assessed the stability and conformational changes of the ERK1-compound complexes over time. Results indicated that Silandrin and Hydroxytuberosone maintained structural similarity to free ERK1, suggesting their stability as potential inhibitors. The study primarily uses in silico techniques, such as molecular docking and molecular dynamics simulations, to determine the strength and binding affinity of the selected inhibitors. These computational methods can only partially replicate the complexity of biological systems, though they offer insightful information. The effectiveness of these compounds in a biological environment will be validated through in vitro or in vivo research. Further validation and exploration of the selected compounds, Silandrin and Hydroxytuberosone, could offer new avenues for AD therapeutics.
Supplemental Material
sj-docx-1-alz-10.1177_13872877241309592 - Supplemental material for Discovering potential ERK1 inhibitors from natural products for therapeutic targeting of Alzheimer's disease
Supplemental material, sj-docx-1-alz-10.1177_13872877241309592 for Discovering potential ERK1 inhibitors from natural products for therapeutic targeting of Alzheimer's disease by Mohammad Taufeeq, Arunabh Choudhury, Afzal Hussain, Mohamed F Alajmi, Taj Mohammad, Anas Shamsi and Md Imtaiyaz Hassan in Journal of Alzheimer's Disease
Footnotes
Acknowledgments
MIH acknowledges the Council of Scientific and Industrial Research for financial support [Project No. 27(0368)/20/EMR-II]. The author acknowledge the generous support from the researchers supporting project number (RSPD2025R980) King Saud University, Riyadh, Saudi Arabia. AS thanks to the Ajman University for providing APC support.
Authors contribution
Mohammad Taufeeq (Conceptualization; Data curation; Investigation; Methodology; Resources; Validation; Writing—original draft); Arunabh Choudhury (Data curation; Formal analysis; Investigation; Resources; Validation; Writing—original draft); Afzal Hussain (Data curation; Formal analysis; Project administration; Resources; Validation; Writing—review & editing); Mohamed F. Alajmi (Data curation; Formal analysis; Funding acquisition; Project administration; Supervision; Writing—review & editing); Taj Mohammad (Conceptualization; Data curation; Formal analysis; Resources; Validation; Visualization; Writing—original draft); Anas Shamsi (Conceptualization; Data curation; Funding acquisition; Methodology; Software; Supervision; Visualization; Writing—original draft); Md. Imtaiyaz Hassan (Conceptualization; Formal analysis; Funding acquisition; Project administration; Software; Writing—review & editing).
Funding
The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work is supported by the Central Council for Research in Unani Medicine (CCRUM), Ministry of AYUSH, Government of India (Grant No. 3-69/2020- CCRUM/Tech). The Research work was financially supported by “Researchers Supporting Project number RSPD2025R980”, King Saud University, Riyadh, Saudi Arabia.
Declaration of conflicting interests
Anas Shamsi and Md. Imtaiyaz Hassan are Editorial Board Members of this journal but were not involved in the peer-review process of this article nor had access to any information regarding its peer-review.
The remaining authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Data availability
The data underlying this article is available within the manuscript.
Supplemental material
Supplemental material for this article is available online.
References
Supplementary Material
Please find the following supplemental material available below.
For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.
For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.
