In Silico Study of the Effect of Thymoquinone on Three Pre-apoptotic Factors of Bad, Bak, and Bim

Background: Apoptosis is one of the most vital mechanisms in the removal of old and damaged cells, especially cancer cells. Many polyphenolic and antioxidant compounds can inhibit the growth and proliferation of cancer cells by inducing apoptosis. This study intended to evaluate the effect of the thymoquinone (TQ) on the three pre-apoptotic factors Bad, Bak, and Bim in the simulation environment. Methods: The Protein Data Bank (PDB) files of three pre-apoptotic proteins Bad, Bak, and Bim were obtained from PDB database and the molecular structure of TQ from PubChem database. The optimization, simulation, molecular docking, and molecular dynamics (MD) computation were performed using AutoDock, VMD, and GROMACS packages on each one of the proteins in free mode and ligand binding mode. Results: The number and type of hydrogen and hydrophobic bonds at the binding site of TQ and preapoptotic factors Bad, Bak, and Bim were detected at the lowest level of energy. The lowest binding energy level of TQ had the highest tendency to bind to the BAD molecule. After the binding of TQ to proteins, the radius of gyration (Rg) of Bim increased while the Rg of Bad and Bak decreased. However, the secondary structures (Turn, Coil, Helix, and Bend) were greatly affected in the binding of TQ to Bad, Bak, and Bim. Conclusion: The variations of binding energy indicated that TQ can affect the three pre-apoptotic factors Bad, Bak, and Bim. This bond seems to increase their activity by variation in the secondary structure of the Bim specific residues.

inhibiting the process, this study tried to investigate the effects of TQ, a powerful flavonoid compound, on three pre-apoptotic factors such as Bad, Bak, and Bim proteins that play an important role in the apoptotic pathways.

Preparation of PDB Files
We obtained the PDB files of the pre-apoptotic proteins (Bad ID: 1G5J, Bak ID: 1BXL, and Bim ID: 1PQ1) from Protein Data Bank server (http://www.rcsb.org). Then, we obtained the molecular structure of TQ (CID: 10281) from PubChem server, converted it to PDB files and optimized by ArgusLab software.

Simulation and Molecular Dynamics Studies
Studies on the molecular dynamics (MD) simulation of the three above-mentioned protein structures were initially carried out in pure water to allow the corresponding structures to be in equilibrium under conditions of constant temperature, pressure, and concentration. The simulations of three apoptotic complexes (Bad, Bak, and Bim) were performed using the GROMACS (version 4.6.1) software package with G43A1 force field in the water solvent. A 140 mM concentration system was prepared by adding the calculated Na and Cl ions. In this study, the SPC216 model was used (18). The output PDB file was then used as an input structure to simulate complexes in molecular docking.

Molecular Docking
Molecular docking studies were done by AutoDock software in Linux operating system with a 64-bit system and Intel (R) Core (TM) i7 CPU Server in Clinical Biochemistry Research Center of Basic Health Sciences Institute of Shahrekord University of Medical Sciences. Molecular docking of TQ on three pre apoptotic factors Bad, Bak, and Bim was done to find the best binding sites for the ligand-receptors and to determine the most stable free energy state of ligand-receptors. In this study, we built a Grid Box with 60×120×72 nm (x×y×z) for proteins, after the production of PDBQ and PDBQT, TQ file as a ligand and Bad, Bak, Bim as receptors. After 200 stages of molecular docking run on the ligand, we used the Genetic Algorithm and Lamarckian Genetic Algorithm parameters. The data obtained from the out put (n.dlg) file were analyzed (3).
In this study, we used LigPlot+ software to specify the number of hydrophobic and hydrogen bonds between TQ and three pre-apoptotic proteins (Bad, Bak, and Bim). Then, the type and number of amino acids presented in the binding site were identified.

Molecular Dynamics Simulation
At the last stage of the simulation of MD, the complex of TQ and all three proteins (Bad, Bak, and Bim) was formed in 140 mM water and salt in accordance with the abovementioned method, and the paths stored in the simulation were used to analyze the structural parameters of the interaction complex. In the following, we used Grapher 10 to analyze the results of the MD simulation of preapoptotic molecules without TQ compared to the results of the MD simulation of the ligand complex with each of the pre-apoptotic molecules (19). The temperature was set to 300 K for all the 10 nanoseconds (ns) of simulation time.

Statistical Analysis
The data were analyzed using SPSS version 22.0 (Chicago, IL, USA). The paired-samples t test was performed to analyze MD. P < 0.001 was considered the significance level.

Results
The results of molecular docking between TQ and three pre-apoptotic proteins Bad, Bak, and Bim are shown in Table 1. The lowest estimated free energy of binding (BE) (-6.04 kcal/mol) and the least estimated inhibition constant (EIC) (37.4 µM) belong to the binding of TQ with  The Effect of Thymoquinone on Pre-apoptotic Factors  Bad molecule. Figure 1 and Table 1 show the hydrogen and hydrophobic bonds in the binding sites of TQ and the preapoptotic proteins (Bad, Bak, and Bim). The hydrogen and hydrophobic bonds between bonded residues and TQ in binding site are shown in Figure 1.
The binding sites of TQ with three pre-apoptotic factors (Bad, Bak, and Bim) are shown in Figure 2.
Root mean-square deviations (RMSD) for Bad, Bak, and Bim in complex with TQ at 10 ns of simulation are shown in Figure 3. Figure 4 shows the total energy (TE) for Bad, Bak, Bim in complex with TQ at 10 ns of simulation. Figure 5 shows the radius of gyration (Rg) for Bad, Bak, and Bim in complex with TQ at 10 ns of simulations. The variations in the secondary structures of Bad, Bak, Bim and the values of these variations as a single molecule and in complex with TQ are displayed in Table 2.

Discussion
The results of this study indicated that TQ with a high binding affinity was able to affect the structure and function of the pre-apoptotic factors Bad, Bak, and Bim. TQ at the binding site provides the most hydrophobic bonds to the Bad molecule, while no hydrogen bonds occurred at the binding site. However, it seems that the Bak protein with two hydrogen bonds and six hydrophobic bonds had more tendency to bind to TQ ( Table 1).
The simulation and MD studies on some antioxidants showed that Gallic acid affected the Bad molecule more compared to the other two molecules, Bak and Bim (19). A study on the effect of carvacrol on the apoptotic factors of Bak, Bax, Bim, Bid, Apaf1, and P38 in a simulated environment showed that carvacrol has been able to activate apoptotic factors in different modes (20). A study on the effects of quercetin on pre-apoptotic factors (Bad, Bak, and Bim) indicated that the number of hydrogenbonds created between quercetin molecules with Bad factor was higher compared to the two Bak and Bim factors. The quercetin docked to Bad molecule had the lowest energy level compared with two other factors. The mean RMSD increased during 10 ns of simulation time for Bad and Bak proteins in the presence of quercetin and decreased for the Bim molecule and it was indicated that Bim had the highest flexibility in the presence of Quercetin compared to free molecule (21).
TQ is a powerful compound that can induce apoptosis  Statistical analysis was done by independent samples t test. Each point represents mean ± SD. The difference between the red and blue lines was significant (P < 0.001). and prevent the growth of cancerous cell (22). The effect of TQ and its antioxidant activity and molecular mechanism on cancer cells have always attracted the attention of researchers (23). Studies have shown that TQ affects the mitochondrial membrane and then releases cytochrome C, which inhibits the activity of the anti-apoptotic factors Mcl-1 and Bcl-xl, and also increases the activity of the apoptotic factors of Bax and AIF (24).
The MD results of this study indicated that all three preapoptotic proteins Bad, Bak, and Bim have been stable at the end of 10 ns in both the pre-docking phase and postdocking phase in TQ binding (Figure 3). Simultaneously with the binding of TQ and activation of apoptotic factors, high amounts of energy have been released by the end of the 10 nanoseconds (Figure 4). The radius of gyration (Rg) of proteins, which is a scale to measure the availability of active site of the protein, is changed by the bonding of TQ to all the three Bad, Bak, and Bim factors. The amount of Rg in binding with the Bim molecule increased, while it decreased with the binding of TQ to the Bad and Bak molecules ( Figure 5). The increase in the amount of Rg in the Bim molecule after binding with TQ seems to have provided conditions for greater activity of this molecule. The results of the MD of this study showed that the binding of TQ to the proteins Bad, Bak and Bim produces changes in the secondary structure of these proteins, which would be very effective in their function ( Table 2).

Conclusion
This simulated and molecular dynamic study shows that TQ with high tendency to bind the three pre-apoptotic factors of Bad, Bak, and Bim can induce the intense changes in their secondary and tertiary structures. These fundamental changes can affect their activity. It is thought that an increase in the amount of radius of gyration causes an increase in protein activity. Apoptotic properties of TQ may be due to its direct effects on apoptotic factors.