Abstract

Objective. MAPK3 activates several nuclear transcription factors, including c-Jun and c-fos, by phosphorylating its downstream cytoplasmic protein, thereby contributing to cell proliferation and survival. Different carcinomas’ initiation, progression, cancer cell metastasis, and drug resistance have been associated with MAPK3 overexpression. Given the need for new and potent MAPK3 inhibitors, this study aimed to explore the potential of anthraquinones (AQs) as organic compounds capable of inhibiting MAPK3. Methods. Using AutoDock 4.0 software, the binding affinity of 21 AQs to the receptor’s active site was evaluated. AQs were ranked based on their ΔGbinding values to the receptor’s active site, with the highest rankings receiving the most favorable scores. The Discovery Studio Visualizer tool was used to demonstrate the interaction modes between the highest-ranked AQs and the MAPK3 catalytic site. Furthermore, a 100-nanosecond molecular dynamics (MD) computer simulation was performed to assess the stability of the docked pose of the most potent enzyme inhibitor identified in this study. Results. The binding affinity of emodin-8-glucoside, aloe-emodin 8-glucoside, pulmatin, rhodoptilometrin, and hypericin to the receptor’s ATP binding cleft was noteworthy, as the ΔGbinding values were < 10 kcal/mol. In addition, emodin-8-glucoside, aloe-emodin 8-glucoside, and pulmatin were found to have inhibition constant values at the picomolar concentration. According to our computer simulation results, the docked pose of emodin-8-glucoside within the active site of MAPK3 achieved a stable state after 70 ns. In other words, the root mean square deviation (RMSD) graph indicated stability within the 70–100 ns timeframe. Conclusion. Inhibition of MAPK3 by emodin-8-glucoside, aloe-emodin 8-glucoside, pulmatin, rhodoptilometrin, and hypericin may have therapeutic potential in cancer treatment.

1. Introduction

Mitogen-activated protein kinases (MAPKs), a family of serine/threonine protein kinases, are highly conserved and have been shown to regulate a variety of cellular processes, including cellular differentiation, proliferation, survival, and apoptosis [1, 2]. It consists of multiple crucial components and phosphorylation events that significantly contribute to tumorigenesis. The activated kinases transmit signals from the extracellular environment, leading to cell growth, proliferation, differentiation, migration, and apoptosis [3]. As a critical signaling hub, the MAPK pathway integrates extracellular signals to regulate cellular differentiation, proliferation, and survival, as well as drug resistance [4, 5]. The pathway is composed of a range of kinases, including MAPKs (MAPK1 (ERK2) and MAPK3 (ERK1)), MEKs, RASs, RAFs, adaptor molecules, and specific negative regulators of ERK1/2 (DUSP3/5/6/7/9). The activation of MAPK1/3, KRAS, HRAS, and BRAF is widely recognized to contribute to several human cancers through transcriptional activation, interactions with other cancer-related pathways (e.g., JAK/STAT and PI3K pathways), and activities that modulate the immune response in some cancers [610]. Around 20% of cases of head and neck squamous cell carcinoma (HNSCC) exhibit activating mutations in the MAPK pathway [11]. Bayat et al. [12] reported that the MAPK1/MAPK3 signaling pathway has a noteworthy involvement in oral squamous cell carcinoma patients who have a poor prognosis, as opposed to those with a favorable prognosis. The interaction between MAPK1 and MAPK3 with a complex network of structures is well established, and their role in promoting the malignant behavior of cancer cells by changing the metabolic signaling pathway is significant [13]. Several reports have defined MAPK3 as an oncogene factor in several malignancies, including breast, ovarian, colorectal, liver, lung, thyroid, and gastric cancers [1418].

The ATP binding site of MAPK3, the kinase active site, is the target of type-I inhibitors. In contrast, type-II inhibitors bind to the allosteric pocket of enzymes, and these two groups make up the primary categories of kinase inhibitors [19].

Anthraquinones (compounds with a 9, 10-anthracene skeleton) are natural substances commonly employed in traditional Chinese medicine [20]. The biological properties of these compounds are multifaceted, encompassing antitumor, analgesic, antibacterial, antimalarial, antioxidant, and anti-inflammatory activities [2123]. The outstanding antitumor activity demonstrated by anthraquinone derivatives has generated considerable interest recently, with inevitable results being authorized for clinical use as antitumor medications. Clinical compounds, including mitoxantrone, doxorubicin, daunorubicin, epirubicin, and emodin, are employed to treat a diverse range of cancers [2427]. In this regard, the antitumor effect of mitoxantrone and doxorubicin has been attributed to their ability to induce apoptosis [28, 29]. The extraction of anthraquinones from Damnacanthus and Morinda spp. represents a promising avenue for obtaining natural anticancer compounds [30].

The present study proposes that AQs serve as effective inhibitors of MAPK3 activity. This inhibition can downregulate downstream signaling pathways, ultimately reducing cell differentiation and survival rate. Hence, a molecular docking analysis was conducted to assess the binding affinity of multiple AQs to the active site of MAPK3. The identification and introduction of the top-ranked MAPK3 inhibitors were based on the ΔGbinding values calculated between the AQs and the enzyme’s active site. The subsequent step involved analyzing the interactions between the top-ranked AQs and residues present inside the MAPK3 catalytic cleft. Subsequently, the stability of the docked pose of the most effective MAPK3 inhibitor was compared with that of a standard drug using molecular dynamics (MD) simulation. The present outcomes hold promise for the development of cancer treatments.

2. Materials and Methods

2.1. Preparation of the Receptor and Ligands for Structural Analysis

Obtaining the 3D coordinates of MAPK3 (PDB ID, 4QTB; X-ray resolution, 1.4 Å) involved downloading from the RCSB database (https://www.rcsb.org, [31, 32]) and visualizing using BIOVIA Discovery Studio Visualizer (DSV) version 19.1.0.18287. Two polypeptide chains (A and B) with the same residues (348; 27–374) were present in the PDB file. Chain A was picked for in silico analyses, and the Notepad++ software was used to remove the water molecules and 38Z from the PDB file. The active site-interacting residues in 4QTB were determined by [1] examining the docked pose of 38Z using DSV software in a two-dimensional view and [2] referencing the publication by Chaikuad et al. [33]. The 4QTB file was subjected to energy minimization to attain its most stable 3D conformation using Swiss-PdbViewer version 4.1.0, which can be accessed at https://spdbv.unil.ch. By employing the GROMOS 43B1 force field, Swiss-pdbViewer can analyze the structure’s energy and correct any irregular geometry via energy minimization. The selection of 21 AQs was made to explore potential MAPK3 inhibitors. An analysis of the MAPK3 ATP-binding cleft involved comparing the binding affinities of the studied AQs with that of ulixertinib (DB13930), a standard drug obtained from the DrugBank database (https://go.drugbank.com/, [34]). As part of our previous investigation [35], we followed an energy minimization procedure on the AQs. By adjusting the MAPK3 3D structure to include Kollman charges and polar hydrogen bonds and applying local charge and rotational motion to ligands structures, PDBQT files were generated for both the receptor and ligands using MGL tools [36].

2.2. Procedure for Docking Analysis

The setup used to conduct the docking analyses was a Windows-based PC with the subsequent specifications: 64-bit system type, 32 GB of installed RAM, and an Intel Core i7 processor [37]. Implementing the AutoDock version 4.0 tool and a semiflexible docking algorithm [38], the ΔGbinding values between AQs, a standard drug, and MAPK3 ATP-binding cleft were computed in kcal/mol units. Our earlier report [19] highlighted the presence of 29 residues in the ATP-binding cleft of MAPK3. As a result, the grid box was set up with the subsequent values for docking analysis: X-dimension, 84; Y-dimension, 60; Z-dimension, 70; X-center, 33.335 Å; Y-center, 55.015 Å; Z-center, 49.3 Å, and spacing, 0.375 Å. To assess the binding affinity between the studied AQs, the positive control inhibitor, and the MAPK3 active site, we constructed several independent docked models for each ligand using the Lamarckian genetic algorithm (GA) approach. The Lamarckian GA docking calculations utilized the following parameters: 50 GA runs, a population size of 150, a maximum number of evaluations set at 2500000, a maximum of 27000 generations, and only the top 1 individual automatically survives. The gene mutation rate was set at 0.02, while the crossover rate was set at 0.8, with a two-point crossover mode. The mean of the Cauchy distribution for gene mutation was set at 0.0, and the variance was set at 1.0. In addition, the worst individual was picked after 10 generations. The docking results were then grouped based on the root mean square (RMS) tolerance of 2.0 Å [36], and the most negative ΔGbinding value within the most significant cluster was chosen for evaluation. Finally, the generation of molecular visualizations was carried out using the DSV software.

2.3. Molecular Dynamics Analysis

The receptor was subjected to MD investigations alone and in complex with [1] the most highly ranked AQ from the molecular docking results and [2] the standard drug ulixertinib. Computer simulations lasting 100 nanoseconds (ns) were employed to conduct MD analyses using Discovery Studio Client software version 16.1.0.15350. A more powerful computer configuration was employed for MD simulations than for docking analyses (system type, 64-bit; installed RAM, 64 GB DDR5; and processor, Intel 24-Core i9-13900KF). The MD simulation was conducted for 100 nanoseconds (100 ns) and implemented using advanced settings. The solvation model utilized was the explicit periodic boundary, while the cell shape was orthorhombic with a minimum distance set at 10 Å from the boundary. The solvent used was water, with a target temperature of 310 K. The force field utilized CHARMm, and the charge distribution was point-based [39, 40]. The root mean square deviation (RMSD) and root mean square fluctuation (RMSF) of the MAPK3 backbone atoms were investigated throughout MD simulations. In addition, the radius of gyration (ROG) and total energy of the receptor were computed to obtain more accurate outcomes.

3. Results

3.1. Binding Affinity Assessment

The top-ranking MAPK3 inhibitors were identified in this study by determining the ΔGbinding value of five AQs with the ATP-binding cleft, which was below 10 kcal/mol. The binding energies for emodin-8-glucoside, aloe-emodin 8-glucoside, pulmatin (chrysophanol-8-0-glucoside), rhodoptilometrin, and hypericin were calculated to be 12.78, 12.42, 12.41, 10.75, and 10.13 kcal/mol, respectively, highlighting their potential as effective inhibitors of MAPK3. The highly potent inhibitory effect of emodin-8-glucoside, aloe-emodin 8-glucoside, and pulmatin on MAPK3 was evident from the estimated Ki value, within the picomolar (pM) range. On the other hand, rhodoptilometrin and hypericin exhibited inhibition at the nanomolar (nM) level, indicating their relatively weaker inhibitory potential. Based on the AutoDock 4.0 tool’s calculation of the ΔGbinding value between ulixertinib and MAPK3, a total of 10 AQs, including top-ranked compounds as well as emodic acid, aloe-emodin, sennidin B, purpurin, and chrysophanol, exhibited a superior binding affinity to the MAPK3 catalytic site than the reference drug. The ΔGbinding value between ulixertinib and MAPK3 was found to be 9.12 kcal/mol (Figure 1). To provide a comprehensive overview of the ligand-MAPK3 affinity, Table 1 displays the calculated ΔGbinding and Ki values between MAPK3 and all ligands. Furthermore, Table 2 presents a detailed analysis of all energy types computed between the ATP-binding site of MAPK3 and the top-ranked AQs. The Gibbs free energy of binding between a ligand and receptor is influenced by several energy factors, including intermolecular energy, internal energy, torsional free energy, and the unbound system’s energy, which have been defined in previous reports [41, 42].

3.2. Evaluating Modes of Interaction

Through the DSV tool, the interactions between the top-ranked AQs, a reference drug, and residues positioned within the MAPK3 ATP cleft were revealed. Computer simulations lasting 100 ns were also conducted to explore the interactions between emodin-8-glucoside, ulixertinib, and residues located within the protein’s active site. Among the ligands studied, emodin-8-glucoside, aloe-emodin 8-glucoside, and rhodoptilometrin formed the highest number of hydrogen bonds (n = 5), whereas rhodoptilometrin displayed the most hydrophobic interactions (n = 6) with the residues of MAPK3. The reference drug formed four hydrogen bonds and one hydrophobic interaction. Before the MD simulation, emodin-8-glucoside, the most effective MAPK3 inhibitor in this research, interacted with the residues of the MAPK3 ATP-binding cleft through five hydrogen bonds and five hydrophobic interactions. Nevertheless, this compound formed seven hydrogen bonds and three hydrophobic interactions after the MD simulation. The interactions between the ligands and the residues of the MAPK3 ATP-binding cleft are presented in Figure 2 and summarized in Table 3. Notably, Table 3 does not include any hydrogen bonds with a distance exceeding 5 Å.

3.3. Robustness of the Docked Conformations

The RMSD range for MAPK3 backbone atoms bound to emodin-8-glucoside was computed between 1.6 and 2.28 Å following 100 ns computer simulation. The RMSD value for a free receptor, however, was found to be between 1.55 and 2.57 Å. In the case of the enzyme complexed with ulixertinib, the minimum and maximum RMSD values were found to be 1.25 and 2.02 Å, respectively. The receptor exhibited more excellent stability when in complex with emodin-8-glucoside in comparison to being free, as suggested by the results. The simulation results showed that the protein’s backbone atoms remained stable after approximately 70 ns of simulation when hindered by emodin-8-glucoside, as indicated in Figure 3(a).

As per the RMSF plots, a notable difference in fluctuation was observed around the 64–124 region. Specifically, the fluctuation of protein complexed with emodin-8-glucoside was lower than that of the free protein (Figure 3(b)). Given the proximity of this region to the MAPK3 active site [19], emodin-8-glucoside is believed to have stabilized the protein’s active site. However, the MAPK3 active site appears to have been more stabilized by the reference drug than by emodin-8-glucoside.

The MAPK3-emodin-8-glucoside total energy was found to be lower than that of free protein during the 100 ns MD simulations. As well, the total energy exhibited by MAPK3-ulixertinib was the lowest compared to that of free MAPK3 and MAPK3-emodin-8-glucoside during the 100 ns MD simulation, according to Figure 3(c).

The results revealed that the ROG value of free MAPK3, MAPK3-ulixertinib, and MAPK3-emodin-8-glucoside increased after the 100 ns MD simulation. As per the data presented in Figure 3(d), the ROG for the enzyme complexed with emodin-8-glucoside was lower than that of free MAPK3 during the first 50 ns computer simulation. Throughout the 50–100 ns MD simulation range, the ROG values of free MAPK3, MAPK3-ulixertinib, and MAPK3-emodin-8-glucoside were in close proximity to each other. The superimposed structures of free MAPK3, MAPK3-emodin-8-glucoside, and MAPK3-ulixertinib, before and after MD simulations, are presented in Figure 4, using the DSV tool. In addition, Figure 5 depicts the incorporation of emodin-8-glucoside in the MAPK3 ATP-binding grove, both pre- and post-MD simulation, carried out with the help of Chimera version 1.8.1.

4. Discussion

The MAPK3 gene encodes an upstream regulator of the MAPK cascade that plays a crucial role in various biological processes associated with apoptosis and cell survival. Aberrant expression of MAPK3 has been linked to the initiation, progression, metastasis, and resistance to the treatment of various human cancers, highlighting the urgent need to explore and develop innovative and potent MAPK3 inhibitors [19]. The objective of this study was to identify natural AQs that could serve as promising inhibitors of MAPK3. The results indicated that the calculated binding free energy between MAPK3 and five AQs (emodin-8-glucoside, aloe emodin 8-glucoside, pulmatin, rhodoptilometrin, and hypericin) was <10 kcal/mol, making them the most potent MAPK3 inhibitors in this study. The recorded ΔGbinding value of 9.12 kcal/mol between MAPK3 and ulixertinib indicates that the top-ranked AQs in this study had a stronger binding affinity to the MAPK3 ATP-binding cleft compared to the standard drug. Furthermore, emodin-8-glucoside, aloe emodin 8-glucoside, and pulmatin exhibited Ki values in the picomolar range.

The MAPK3 active site exhibited a ΔGbinding value of 12.78 kcal/mol upon binding with emodin-8-glucoside. Furthermore, the binding energies between aloe emodin 8-glucoside, aloe-emodin, emodin, and MAPK3 were estimated to be 12.42, 9.63, and 8.85 kcal/mol, respectively. These results suggest a sugar moiety enhances emodin’s binding affinity to MAPK3. Prior to undergoing MD simulations, emodin-8-glucoside was shown to establish five hydrogen and five hydrophobic interactions with residues Tyr53, Val56, Lys71, Glu88, Gln122, Asp123, Cys183, and Asp184 situated within the ATP-binding cleft of MAPK3. This compound exhibited seven H-bonds and three hydrophobic interactions with Glu50, Val56, Lys71, Met125, Lys131, Ser170, Leu173, and Cys183, after MD simulations.

The use of Chinese herbs to treat different cancers has been established for a long time due to their confirmed benefits and negligible side effects. Emodin (1,3,8-trihydroxy-6-methylanthraquinone) is obtained from the root and rhizome of Rheum palmatum L., a medicinal herb with a long-standing tradition [43]. Previous research has demonstrated that emodin possesses inhibitory potential against several cancers, such as hepatocellular carcinoma [44], breast cancer [45], cervical cancer [46], ovarian cancer [47], and bladder cancer [48]. According to Manimaran et al. [49], emodin administration at a dose of 50 mg/kg b.w. was effective in inhibiting the activation of Akt, ERK, P38 MAPK, and DNA methyl transferase (DNMT) in golden Syrian hamsters with 7,12-dimethylbenz[a]anthracene (DMBA)-induced oral carcinoma. This was evident from the downregulation of these markers as observed through Western blotting. Furthermore, emodin was found to alleviate the severity of precancerous lesions, such as dysplasia, in DMBA-treated hamsters. As well, the findings of Lin et al. [50] indicate that aloe emodin can suppress the expression of matrix metalloproteinase-2 via the P38 MAPK-NF-kappa B signaling pathway, thereby inhibiting the invasion of nasopharyngeal carcinoma cells.

According to the findings, the binding of chrysophanol-8-0-glucoside (pulmatin) to the MAPK3 catalytic domain resulted in a Ki value of 793.77 pM and a binding energy of −12.41 kcal/mol. The MAPK3 showed less binding affinity towards chrysophanol than its glycosylated form, with a recorded ΔGbinding value of 9.19 kcal/mol and a Ki value of 184.75 nM. The MAPK3 active site was found to interact with chrysophanol-8-0-glucoside through four hydrogen and four hydrophobic interactions with specific residues, including Tyr53, Val56, Lys71, Gln122, Asp123, and Asp184, before MD simulations.

Chrysophanol, an anthraquinone metabolite obtained from the Rheum genus, has been found to possess anticancer properties in recent studies [5153]. In addition, it exhibits anti-inflammatory activity [54] and provides neuroprotection effects [55]. Rheum genus contains a higher amount of chrysophanol-8-O-glucoside, which is a glycosylated form of chrysophanol, compared to chrysophanol [56]. Our bioinformatics results are in line with those of the research of Kwon et al. [57], who have found that the Cassia tora seed extract, a rich source of chrysophanol, can deactivate the JNK/P38 MAPK signaling pathway, thereby suppressing heat-induced lipogenesis in human sebocytes.

Upon conducting a comparison of the “emodin-8-glucoside”-MAPK3 contacts between the initial docking and final MD-simulated results, several noteworthy points emerged. These are as follows:(a)The hydrophobic interaction between Val56, Cys183, and emodin-8-glucoside remained stable before and after the MD simulation.(b)Cys183 formed a hydrogen bond with emodin-8-glucoside in addition to the hydrophobic interaction after the MD simulation.(c)The contact between emodin-8-glucoside and Lys71 was initially a hydrophobic interaction, while after the MD simulation, it was observed as a hydrogen bond.(d)Following the MD simulation, new residues interacted with the emodin-8-glucoside, including Glu50, Met125, Lys131, Ser170, and Leu173.

Furthermore, the results indicated that the interacting residues between MAPK3 and ulixertinib underwent a complete change following the MD simulation. Prior to the simulation, Tyr53, Asp128, Ser170, and Asn171 were found to interact with ulixertinib. However, after the simulation, hydrophobic and hydrogenic interactions were observed between Tyr130, Lys131, and Lys168 and the reference drug. It is worth noting that, based on the present findings, the interactions between emodin-8-glucoside and MAPK3 were comparatively more stable than those observed with the standard drug.

5. Conclusion

Based on the current study, emodin-8-glucoside, aloe-emodin 8-glucoside, pulmatin, rhodoptilometrin, and hypericin possess strong inhibitory properties against MAPK3. Among these metabolites, emodin-8-glucoside, aloe-emodin 8-glucoside, and pulmatin are the most potent, with Ki values at the picomolar scale. The computer simulation indicated that the 3D conformation of MAPK3 remained stable when in complex with emodin-8-glucoside for approximately 65 ns. These findings may offer valuable insights to researchers aiming to create novel drug therapies for numerous types of cancer. However, additional research is needed to confirm the present results, including more extended MD simulations and in vitro and in vivo validation experiments.

Data Availability

The datasets used and/or analyzed during the current study are available from the corresponding author upon reasonable request.

Ethical Approval

The present study was approved by the Ethics Committee of the Hamadan University of Medical Sciences, Hamadan, Iran (ethics no. IR.UMSHA.REC.1401.254).

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Authors’ Contributions

AT and SV-A designed the study. MM conducted docking operations. AT performed MD simulations. AT and SV-A analyzed and discussed the results. AT wrote the manuscript. SV-A edited the manuscript. All the authors read and approved the final version of the manuscript.

Acknowledgments

The authors thank the Research Center for Molecular Medicine, Dental Research Center, Hamadan University of Medical Sciences, Hamadan, Iran, for their support.