The X-ray crystallographic structure of p38α MAP kinase was downloaded from the RCSB Protein Data Bank (http://www.rcsb.org/pdb/home/home.do). The structure was cocrystallized with the pyrimidine-urea analogue 15 (PDB ID 2GHL, 2.10 Å X-ray resolution). The structure was stripped of the pyrimidine-urea analogue, where after removing water molecules, the charge and polar hydrogen atoms were added. The incomplete residues were corrected and the energy minimization of the modified structure was set to the default cut-off root mean square deviation (RMSD) value of 0.3 Å using OPLS 2001 force field by SPDB viewer (27). The structure of curcumin and its derivatives were retrieved from the PubChem database (Fig. 1).

Molecular docking simulation was carried out using AutoDock 4.2.6 tool (28). This tool uses binding free energy evaluation for mining the top compound-protein pose. The tool uses binding site definition to create a rectangular grid. The execution of the docking runs are performed within the grid. The energy values calculated by this suite are the culmination of the intermolecular energy of the protein, the internal energy of the ligand and the torsional free energy. The default parameters were set, and the maximum evaluations were performed to elucidate the top 10 binding modes for each run. The best compound among curcumin and its derivatives was selected for further evaluation at the atomic level using molecular dynamics simulation. The results were analyzed using Discovery Studio Visualizer (29).

GROningen MAchine for Chemical Simulations GROMACS 4.5.6 atomic simulation tool was used to evaluate the ligand-p38α complexes over 100 nanoseconds (nsec) (30). GROMOS 43a1 force field was used for p38α and the force field for the top ligand was calculated using PRODRG webserver (31,32). The complexes were kept in a water bath, where water molecules were represented by simple point charge (SPC216) model. The system was ensured to be overall charge neutral before being subjected to energy minimization process, where position restrain procedure was performed with NVT and NPT ensembles, both of which follow Berendsen coupling scheme. Time duration of 1 nsec with coupling constant of 0.1 psec at a constant temperature of 300 K was used for NVT ensemble. This was followed by NPT ensemble with a time duration of 1 nsec; here a constant pressure of 1 bar was subjected to a coupling constant of 5 psec. The complexes were then subjected to simulation for 100 nsec at 300 K. For each run, the SHAKE algorithm for bond lengths and particle mesh Ewald method for long-range electrostatic interactions was employed, where a 14 Å cut-off for van der Waals interactions, a 12 Å cut-off for Coulomb interaction with updates every 2 steps (30,33). The trajectories generated were subjected to GROMACS inbuild tools for evaluation. The results generated were visualized using Pymol, Ligplus and VMD (34-36) and the graphs were plotted by using Grace program (37).

The trajectories generated were used for MM-PBSA calculations (38,39). For each trajectory 1001 complexes were retrieved, each after 100 psec. Binding modes were analyzed using Grace program (37). The tool gave us the binding free energy of the four complexes. The binding free energy (ΔG_bind) calculated by the tool includes the following species: Δ G bind = G complex − G protein − G ligand = Δ E MM + Δ G sol − T Δ S (I) Δ E MM = Δ E val + Δ E ele + Δ E vdw (II) Δ G sol = Δ G p + Δ G np (III) Δ G np = γ SASA + β (IV)

In the equations G_complex, G_protein and G_ligand represent the free energy of the respective species, ΔE_MM, ΔG_sol and TΔS represent gas phase energy, solvation energy and an entropy term. The gas phase energy is the summation of the internal energy of the bonds, angle and torsion (ΔE_val), electrostatic interaction energy (ΔE_ele), and van der Waals interaction energy (ΔE_vdw). The salvation free energy comprises the polar salvation free energy (ΔG_p), and the nonpolar salvation free energy (ΔG_np).

g_covar, g_anaeig in-build tools of GROMACS suite were used for indentifying the dominant motions of the selected complexes. PCA technique was used to trace these high-amplitude concerted motions in the p38α-curcumin-based complexes (40). Eigenvectors were generated from the covariance matrix generated by simple linear transformation. The first two eigenvectors (EV1 and EV2) having cosine content <0.2 were used to study free energy landscape (FEL). g_sham in-build tool was used to calculate FEL.

THP-1 monocytes were procured from American Type Culture Collection (ATCC; Manassas, VA, USA) and cultured in RPMI-medium (Sigma-Aldrich, St. Louis, MO, USA) supplemented with 10% FBS and 100 U/ml PenStrep (Gibco, Grand Island, NY, USA). The cells were grown under standard culture conditions at 37°C under 5% CO₂ in a humidified incubator. For macrophage differentiation, monocytes were cultured with macrophage CSF (M-CSF) at 100 ng/ml. Stock solution of 100 mM cyclocurcumin (Sigma-Aldrich) was prepared in RPMI-medium and diluted to different concentrations. Cells were stimulated with 10 ng/ml lipopoly-saccharide (LPS) (Sigma-Aldrich) for TNF-α production.

3-(4,5-Dimethylthiazol-2-yl)-2,5-diphenyltetrazolium bromide (MTT) assay is a standard test to measure cell viability and is based on the conversion of MTT to formazan crystals by active mitochondrial dehydrogenases. Briefly, a cell suspension containing ~3×10⁴ cells was plated into each well of a 96-well plate and allowed to attach for 24 h. Cyclocurcumin prepared in DMEM was added to the wells at a concentration of 10, 20, 40, and 60, 80 and 100 µM for 24 h. After treatment, 100 µl of MTT (Sigma-Aldrich) solution was added to each well at a concentration of 0.1 mg/ml (dissolved in PBS). After 4 h of incubation at 37°C in dark, the solution was removed and 100 µl of dimethyl sulfoxide (DMSO) was added to solubilize the produced formazan. Absorbance was measured at 570 nm using an automated microplate reader (BioTek Instruments, Winooski, VT, USA). The results are expressed as percentages relative to the controls. The percentage of proliferation inhibition was calculated as = (1 - OD_sample/O_Dcontrol) ×100%.

Following treatment with cyclocurcumin and stimulation of human macrophages with LPS, the supernatants were harvested. Concentrations of TNF-α (BD Pharmingen, San Diego, CA, USA) were determined by ELISA, following the manufacturer's instructions. Absorbance was read at 450 nm on an ELISA plate reader (BioTek Instruments) using an in-built software program.

Statistical analysis was carried out by one-way ANOVA and the level of significance was tested at a P-value ranging from 0.01 to 0.001. Data are expressed as the mean ± SD of three independent experiments.

Molecular docking analysis was carried out using AutoDock 4.2 tool; 100 conformations for each compound were generated and the binding energy (kcal/mol) was set as the criteria for the cut-off for the molecular dynamics simulation analysis. Table I shows the results generated by docking curcumin and its selected derivatives against the active site of p38α. The top 4 compounds with the least binding energy were curcumin, bisdemethoxycurcumin, cyclocurcumin and demethoxycurcumin. Protein/ligand complexes from the molecular docking were subjected to post docking analysis, where weak interactions involved in the complex were studied (Fig. 2). Four complexes shortlisted for this analysis showed the favorable binding energy with a range -5 to -7 kcal/mol (Table II). The effective interactions studied here were supported by energies including final intermolecular energy, electrostatic energy and van der Waals energy. The effective interactions of the curcumin-p38α complex were evaluated using Discovery Studio Visualizer 4.5. The binding pocket of the complex involved GLY110, ASN115, ASP112, ASP168, VAL38, LEU104, GLU71, LYS53, VAL105, LEU75, LEU86, VAL30, THR106 and ALA111. The binding energy of curcumin-p38α was calculated to be the highest among the studied complexes, with no non-covalent interaction observed (Fig. 3A). The binding energy of the bisdemethoxycurcumin-p38α complex at -6.38 kcal/mol was the second highest calculated among the studied complexes. The bisdemethoxycurcumin-p38α complex formed one hydrogen bond (Fig. 4A). The interacting residue PHE169 of p38α is involved in active p38α inhibition. The oxygen atom of PHE169 formed a physical non-covalent interaction with bisdemethoxycurcumin, and the distance of the hydrogen bond formed was 2.36 A°. The binding pocket (Fig. 3B) of the bisdemethoxycurcumin-p38α complex involved LEU75, LEU74, THR106, ASP168, GLY170, GLU71, ILE84, PHE169, LEU167, MET109, GLY110, ASP112, ALA111 and ALA157. The third complex was that of cyclocurcumin with p38α protein kinase. The binding energy of the complex was calculated to be −6.12 kcal/mol. In the cyclocurcumin-p38α complex the results showed two hydrogen bonds (Fig. 3C) between the protein and ligand atoms. The interaction residues included LEU104 and VAL105 involved in active p38α inhibition. The binding pose of hit CID 69879809 clearly suggested effective inhibition. The atom 'O₂' of the key residue LEU104 interacted with cyclocurcumin through hydrogen bond interactions (Fig. 4B). The key residue VAL105 of the active binding site cavity also interacted with cyclocurcumin oxygen atoms. Overall, the number of hydrogen bonds confirmed the high affinity of cyclocurcumin towards inhibition of p38α. In the case of the demethoxycurcumin-p38α complex, the results showed 3 hydrogen bonds between p38α and demethoxycurcumin (Fig. 3D). The interaction residues included SER154 and ASN155 (Fig. 4C). The interaction of an 'O' atom of SER154 and two 'O' atoms in the demethoxycurcumin aromatic ring showed a strong interaction pattern. The 'O' atom in both sides of the demethoxycurcumin aromatic ring was more potent in p38α inhibition.

After identification of the top 4 curcumin-based compounds from the molecular docking, molecular dynamic simulation (MDS) was employed to refine the top 4 complexes. MDS is a powerful technique to monitor protein-ligand complexes over a period of time. MDS run at a 100-nsec production was carried out under GROMOS force field. The analysis was evaluated by group properties, RMSD, hydrogen bond interactions and advance properties including PCA and MM/PBSA free energy calculations using trajectory files. All the 4 complexes were stable throughout the run. For each complex, RMSD was evaluated for the convergence of the protein structure towards an equilibrium state after lead binding. Fig. 5 shows the RMSD plot based on backbone atoms. Curcumin and bisdemethoxycurcumin in complex with p38α exhibited a jump at 30 nsec and equilibrated well and both of the complexes converged with an RMSD range almost >0.2 nm. The RMSD value of the bisdemethoxycurcumin-p38α complex was least among the complexes. With an RMSD of 0.3 nm, cyclocurcumin showed the highest deviation. From the RMSD plot of the demethoxycurcumin-p38α complex, the structure equilibrated at 30 nsec and then started to converge near 0.25 nm. The change in RMSD in all 4 complexes up to 30 nsec was due to ligand binding with the active site and after 30 nsec all of the structure obtained stability with less structural deviation at the end of the 100-nsec simulation. After ligand binding, the structural deviation of the bisdemethoxycurcumin-p38α complex was the most stable among the 4 under consideration. Inter-hydrogen bond interactions between protein and ligands were evaluated for the p38α-curcumin-based complexes (Fig. 6). In the case of the curcumin complex (Fig. 6A), NH plot results showed a range of 4–6 hydrogen bond interactions observed throughout the 100-nsec simulation and a maximum of 8 hydrogen bonds. NH analysis confirmed strong inhibition of p38α by curcumin in dynamic system similar to the docking results inferred with 6 hydrogen bonds. In the case of the bisdemethoxycurcumin complex (Fig. 6B), the results showed a range of 3–4 hydrogen bond interaction found throughout the 100-nsec simulation and a maximum of 6 hydrogen bonds similar to the docking results. In the case of the cyclocurcumin complex (Fig. 6C), the results showed a range of 1–5 hydrogen bond interactions found throughout the 100-nsec simulation and a maximum of 5 hydrogen bonds similar to the docking results. The inter-hydrogen bond interaction pattern suggested the plausible mode of the strong binding of the lead candidates with p38α which favored the inhibition mechanism.

The best 4 ligands in complex with p38α from MDS were used to calculate the binding free energy using MM-PBSA method. Snapshots were extracted at every 10-psec of stable intervals from 100-nsec MD trajectory and served as input for calculation. The binding free energy and its corresponding components obtained from the MM/PBSA calculation of the p38α curcumin-based complexes are listed. The results indicated that 4 curcumin-based compounds possessed a negative binding free energy of −141.727, −132.758, −151.455 and −133.56 kcal/mol (Table III). Moreover, van der Waals, electrostatic interactions and non-polar solvation energy negatively contribute to the total interaction energy while only polar solvation energy positively contributes to total free binding energy. The relative binding free energies of the 4 p38α-curcumin-based complexes supported the strong binding in the dynamic system. The interaction with the active site of p-p38α is shown in Fig. 7.

Principal component analysis or essential dynamics was performed on the cyclocurcumin complex with p38α. The analysis was used to monitor the overall strenuous motions of the complex. The results showed that the first 6 eigenvectors formed 90% of the total motion of the complex. The first 2 components of the eigenvectors were projected at 300 K and analysis of this revealed the clusters of stable states for the complex (41). The FEL was plotted over this projection, which provided a clear description of the stability of the complex (Fig. 8).

TNF-α is a key factor in a variety of inflammatory diseases. Keeping in view the role of p38α in regulating TNF-α expression, we evaluated its expression in the culture supernatants of human macrophages (which constitute the major TNF-α-producing cells in highly relevant inflammatory disorders) treated with non-cytotoxic doses of cyclocurcumin. In order to determine the suitable dosage of cyclocurcumin that result in less cytotoxicity, macrophages (obtained with treatment of THP-1 monocytes with CSF), were treated with different concentrations of cyclocurcumin (10–100 µM) for 24 h, and the viability was measured by MTT assay. The results demonstrated that treatment with cyclocurcumin at concentrations such as 10, 20 and 40 µM exhibited the least cytotoxic effects (Fig. 9). Therefore, these concentrations of cyclocurcumin, without cytotoxic effects on macrophages were used. To test the effect of cyclocurcumin, macrophages were treated with the indicated concentrations of cyclocurcumin or 10 µm SB203580, a selective inhibitor of p38 MAPK. This was followed by stimulation with 10 ng/ml LPS for 18 h to induce the production of TNF-α. It was observed that treatment of the LPS-stimulated macrophages with cyclocurcumin led to significant inhibition in the release of TNF-α in a dose-dependent manner, similar to the inhibition caused by SB203580 (Fig. 10). Thus, these data confirm the role of cyclocurcumin in overcoming p38α-induced production of TNF-α and hence can be used as a therapeutic agent to target RA.

In conclusion, natural products have been described as a gold mine for arthritis treatment. These products can suppress the expression of TNF-α, IL-1β, cyclooxygenase-2, lipooxygenase, matrix metalloproteinases or adhesion molecules, or suppress the activation of nuclear factor-κB (NF-κB). All are factors linked to the treatment of RA (42). Curcumin has long been reported to have an antirheumatoid effect. Initial studies have linked it with the downregulation of NF-κB (43). A plethora of research on curcumin has given us an opportunity to investigate it and its derivatives as potential kinase inhibitors for p38α. To the best of our knowledge, our study is the first to use an advanced computational approach to identify the possible DMRDs. MM-PBSA is a potent tool in the field of drug design (44) and this approach was employed in this study for short listing of compounds. From the molecular docking and molecular dynamic simulation studies, the potential binding mode of cyclocurcumin with p38α with stability was revealed as a top compound for the treatment of RA. Finally, inhibition of the release of TNF-α from LPS-stimulated macrophages by cyclocurcumin treatment confirms its role as a potent p38α inhibitor.