Theoretical investigations on the antioxidant potential of a non-phenolic compound thymoquinone: a DFT approach

Thymoquinone (TQ) is a natural compound present in black cumin which possesses potent antioxidant activity without having any phenolic hydroxyl group which is responsible for antioxidant activity. In the present study, computational calculation based on density functional theory (DFT) was executed to assess systematically the antioxidant behavior of this compound by considering geometrical characteristics, highest occupied molecular orbital - lowest unoccupied molecular orbital (HOMO-LUMO), and molecular electrostatic potential (MEP) surface. Thermochemical parameters correlated to the leading antioxidant mechanisms such as hydrogen atom transfer (HAT), single electron transfer-proton transfer (SETPT), and sequential proton loss electron transfer (SPLET) were studied in gas and water media. In addition, the changes of thermochemical parameters such as free energy change (∆G) and enthalpy change (∆H) were computed for hydrogen abstraction (HA) from TQ to hydroxyl radical in gas and water phases to investigate its free radical scavenging potency. The low and comparable values of bond dissociation enthalpy (BDE), proton dissociation enthalpy (PDE), ionization potential (IP), proton affinity (PA), and electron transfer enthalpy (ETE) revealed the antioxidant activity. The ∆G and ∆H also indicated apposite thermodynamic evidence in favor of antiradical capability of TQ. The attack of the free radical occurred preferentially at 3CH position of the molecule.


Introduction
Free radicals form in living cells via biochemical reactions, but some external factors such as radiation, pollutants, physical stress, pesticides, and several treatments also cause its formation [1,2]. The highly reactive free radicals can initiate chain reactions interacting with proteins, lipids, and DNA, and thereby damage the cell tissues [3]. Foods containing antioxidant compounds play a crucial role to check the human illnesses caused by the detrimental effects of free radicals [4]. Polyphenols, phenolic acids, and flavonoids either synthetic or natural can protect biomolecules against unexpected oxidative damage being triggered by free radicals. Thus, they inhibit the oxidative mechanisms that lead to degenerative diseases [5,6].
Thymoquinone (TQ) is an active component of black cumin (Nigella sativa L.) essential oil possessing antineoplastic, anti-inflammatory, neuro-protective, and hepato-protective properties. These properties are attributed to the antioxidant potential of TQ, which seems to be unlikely due to its structure [7]. It has a number of pharmacological properties in which antioxidant activity is significant one [8]. It is well-known that only the compounds having phenolic OH groups are responsible for antioxidant activity and/or free radical scavenging power [9][10][11][12]. Recently, it has been proved that compounds without phenolic OH group can have antioxidant or antiradical capability. Baschieri et al. explained the antioxidant activity of three common nonphenolic terpenoids [13]. Several investigations revealed that some CH bonds contribute to the antioxidant characteristics of natural products [14,15]. The antioxidant potential of TQ is well-recognized in various literature [16][17][18][19], although it does not contain any phenolic OH group.
In the computational method, the antioxidant activity of a compound can be thermodynamically evaluated by some physicochemical descriptors such as bond dissociation enthalpy (BDE), proton dissociation enthalpy (PDE), proton affinity (PA), ionization potential (IP), electron transfer enthalpy (ETE), hardness, softness, and electronegativity. Recently, computational techniques, especially DFT method, have been commendably applied to determine these thermodynamic parameters of antioxidant molecules [3,20].
In the modern computational platform, antioxidant ability of a compound can be successfully determined by monitoring some specific chemical reaction mechanisms. Among them, hydrogen atom transfer (HAT) is familiar one and defined as follows: which relates to the bond dissociation enthalpy (BDE) of the corresponding (A-H) bond that presents direct relation to experimental determination of antioxidant activity [21].
Another important mechanism is the single electron transfer-proton transfer (SETPT) mechanism in which the IPs of the molecule and PDE of AH •+ are accounted [3].
The third pathway is sequential proton loss electron transfer (SPLET) mechanism to investigate the antioxidant potential of a molecule [22]. SPLET is given by: The reactivity and stability of a molecule are considered as significant factors to determine the antioxidant potency which depend upon the difference (ΔE) between the lowest unoccupied molecular orbital (LUMO) and the highest occupied molecular orbital (HOMO) [23]. Large ΔE indicates low reactivity and high stability of the compound. The chemical reactivity also rest on the synergistic effects of different parameters. The E HOMO energy of a molecule shows the electron-donating ability, whereas E LUMO characterizes electron-accepting ability. Ionization potential (IP) and electron affinity (EA) are also roughly accompanying with E HOMO and E LUMO according to Janak's Theorem [24].
Hardness (η) is defined as resisting power toward the polarization of electron cloud of a chemical species [25].
Softness is just a reciprocal of the hardness [26].
Electronegativity (χ) can be described as the ability of an atom to have more attraction to covalently bonded shared electrons toward itself. It is computed as follows: Thus, the antioxidant power of a compound can be determined by taking account all the possible physicochemical characteristics such as geometrical, thermochemical, and orbital characteristics.
However, so far, no or less computational study on the antioxidant power of TQ molecule ( Fig. 1) has been reported until date. The objective of this article was to apply the computational methods to evaluate the antioxidant ability of TQ with its non-phenolic structure. In addition, the stability of the formed radicals from the entitled compound was explained to conduct the comparative radical scavenging ability of the groups undertaken in this study.

Computational methods
All calculations were carried out by using the Gaussian 09 software tool [27]. The structures of TQ, its ions, and radicals were optimized using DFT with the B3LYP method [28][29][30][31] consisting of Becke's three-parameter exchange functional and Lee, Yang, and Parr (LYP) correlation functional, and basis set 6-311+G(d,p) [32] in gas phase and water medium. This methodology demonstrated good agreement between theoretical and experimental outcomes [33,34]. B3LYP was chosen for being one of the most utilized functional tools in computational chemistry, and thus facilitates comparison with upcoming results [35][36][37][38]. All the optimized structures were accepted through vibrational frequency analysis (no imaginary frequency found for optimized structure). The solvation effect of water (ϵ = 78.4), as a universal solvent, on the antioxidant potential of TQ was taken into account employing the conductor-like polarizable continuum model (C-PCM) at the same level of theory. In this model, the solute is positioned in a cavity formed by combination of a series interlocking atomic spheres. Within the cavity site, the dielectric constant is identical as in vacuo; outside, it takes the value of the expected solvent. It assumes that the solvent has an infinite dielectric constant and this method is suitable for highly polar solvents [39]. The thermochemical parameters (BDE, IP, PDE, PA, ETE) associated with antioxidant mechanisms were calculated for the TQ molecule [40,41]. Here, in the equations, TQ has been written as (TQ-H) for better understanding. In order to calculate the reaction enthalpies, it needs the enthalpies of hydrogen atom, proton, and electron. The gas phase values were taken from previous calculations [42,43]. Enthalpies of proton and electron in water were taken from reference [44], and that of hydrogen atom from reference [45].
Changes of Gibbs free energy (ΔG) and enthalpy (ΔH) of the hydroxyl radical scavenging reaction by TQ were calculated according to the simple equations of thermochemistry as follows: Gibbs free energy changes for the reactions are } respectively, and similarly enthalpy changes for the above reactions are Where G(X) and H(X) are the free energy and enthalpy respectively obtained from the optimized geometry of the associated species X calculated at the same level of theory.

Spin density analysis
Spin density distribution is an important quantum property in the estimation of free radical scavenging power of a compound. The stability of a free radical largely depends upon the electron spin density distribution over the radical. Higher delocalization of spin density in the radical reflects its easier formation as well as higher stability [46].

Geometry optimization and conformational analysis
To know the antioxidant behavior of a compound, it is highly recommended to examine its geometrical conformation. For this reason, the geometrical structures of TQ, its radicals, cation, and anions were optimized in the gas and water at B3LYP/6-311+G(d,p) condition using DFT. Some characteristics such as bond length, bond angle, and torsion angle of the entitled compound are tabulated in Table 1. The optimized structure with atom numbering and the schematic projection of TQ are shown in Fig. 1. Very little or no change in the geometrical parameters was observed for optimized TQ molecule in gas phase and water media. Energy changes along the optimization process are depicted in the figure (S- Fig. 1). Optimized structure of TQ in aqueous phase contains 5.39 kcal/mol −1 less energy than that in gas phase.
After withdrawal of proton from the 3CH bond, permanent carbon-carbon double (C=C) bond formed between C3 and C4 with delocalization of π-electrons around C4-C8-C10 of the ring (S- Fig. 3) and in case of 3CH-radical, delocalization of electrons is distributed over C3-C4-C8-C10 areas both in gas and water phase (S- Fig. 2). One important change observed after removal of proton and hydrogen atom from 3CH group was that, the 5CH 3 and 6CH 3 groups came to the same plane of the ring that revealed the whole structure became planer (S- Fig. 2 and S-Fig. 3). Due to the subtraction of proton from 12CH 3 , carbon-carbon double bond was observed between C12 and C9 with delocalization of electrons around C9-C11-C7 of the ring in both phases studied, and similar changes were found for 12CH 3 radical. These delocalization of electrons recommend in favor of antioxidant potential of the studied molecule [41]. But such significant changes were not observed for the elimination of proton and hydrogen atom from 5CH 3 and 6CH 3 groups of TQ.

Bond dissociation enthalpy and spin density: HAT mechanism
Bond dissociation enthalpy (BDE) plays a vital role in the evaluation of antioxidant potential of a compound as it describes the hydrogen atom-donating ability and formation of a stable radical. The lower the BDE, the higher the hydrogen atom-donating tendency as well as the higher free radical scavenging power [47]. The BDEs for concerned (C-H) bonds of TQ were calculated in gas and aqueous phases at B3LYP/6-311+G(d,p) level in DFT and presented in Table 2. From the calculated values of BDE for the possible functional groups of TQ, it was found that 3CH, 12CH 3 , 5CH 3, and 6CH 3 possess BDEs of 77.40, 83.68, 100.25, and 100.19 kcal/ mol −1 in gas phase and 75.80, 83.07, 100.15, and 100.26 kcal/mol −1 in water medium respectively. The values of BDE of the groups in the water medium have been decreased slightly which shows very good agreement with a previous study [3]. The values of BDE of group 12CH 3 and 3CH were less than or close to the literature values computed for other compounds such as ascorbic acid, gallic acid [3], myricetin [48], and myricetin 3,4-di-O-α-L-rhamnopyranoside [49]. The comparable values of BDE of HAT reaction calculated for TQ suggest its antioxidant potential through free radical scavenging mechanism.
In order to elucidate the changes in BDEs, the distributions of spin density of the corresponding TQ • radicals were figured out in the gas phase and water phase, and presented in Fig. 3. It is a reasonably trustworthy quantity to explain the stability of the free radical. The radical with higher delocalization of spin density is more stable [50]. As displayed in Fig. 2, after deduction of a hydrogen atom from the functional groups studied herein, 3CH and 12CH 3 groups had the most delocalized spin density in the corresponding radical which are in good concordance with their sequence of BDE values.
The spin densities for 3CH and 12CH 3 radicals are more distributed, whereas those for 5CH 3 and 6CH 3 are more localized around the concerned groups. This can also reflect the evidence for the more free radical scavenging strength of the TQ molecule with respect to 3CH and 12CH 3 groups. This tendency is almost similar for all the radicals in both gas and water media. As can be concluded from the above discussions, the 3CH and 12CH 3 groups mainly contribute to the Ionization potential and proton dissociation enthalpy: SETPT mechanism SETPT mechanism comprises two stages: donating electron to a free radical R • and subsequently transfer a proton to the anion. In order to justify the possibility of these reactions, it is required to determine ionization potential (IP) for the former step and proton dissociation enthalpy (PDE) for the later one. SETPT mechanism combines IP and PDE to define the antioxidant activity of a compound. The lower IP and PDE values represent higher antioxidant tendency through the SETPT scheme [41]. Table 3 represents the IP of TQ molecule and PDEs of its different CH bonds. The PDE value was the lowest (193.04 and 220.88 kcal/mol −1 in gas and water respectively) for CH bond of 3CH group and the highest at 215.88 kcal/mol −1 of 5CH 3 in gas and 245.34 kcal/mol −1 of 6CH 3 in water. The PDE value of the CH bond of group 3CH was less than that of literature value of other substances demanding antioxidant capacity [51,52]. The lower value of PDE of 3CH bond might be due to the two electrondonating methyl groups attached to the carbon atom. But here, the IP values of TQ were found to be a bit higher in comparison to literature [3] because it was slightly difficult for TQ to donate electron. In fact, the electron-donating ability is related to an extended electronic delocalization over the entire molecule. A molecule having a high degree of π-delocalization is more active [53]. From the above discussion, it can be said that TQ is less likely to pose SETPT formalism of antioxidant behavior.

Proton affinity and electron transfer enthalpy: SPLET mechanism
The third mechanism of the antioxidant activity (SPLET) is characterized by the proton affinity (PA) value which corresponds to the enthalpy of proton dissociation from the neutral molecule, and the electron transfer enthalpy (ETE) related to The values of ascorbic acid and Gallic acid were taken from the reference [3] Gas Water  Table 4. The PA values for the studied CH bond of TQ were comparable to the literature value [3]. PA values of TQ were lower in water than in gas phase which is consistent with previous literature [41].  The values of ascorbic acid and Gallic acid were taken from the reference [3]  The values of ascorbic acid and Gallic acid were taken from the reference [3] The PA values were in the range of 342-381 kcal/mol −1 in gas phase and 295-339 kcal/mol −1 in water. The minimum (342.35 kcal/mol −1 ) PA was found for the CH bond of the 3CH group and maximum (381.73 kcal/mol −1 ) of 6CH 3 group in gas phase. The lower value of PA for 3CH can be due to the electron-pushing character of two methyl groups attached herewith. The ETE values of various CH bonds of TQ were lower in gas than those in water which is also in agreement with literature [53]. The ETEs were in the range of 32-56 kcal/ mol −1 in gas, whereas ETEs were in the range of 74-106 kcal/ mol −1 in water. This also shows very good agreement with the literature of similar study [41]. These values were less than or close to the corresponding values of another similar research [50,53]. The lowest PA values were observed for the 3CH group both in gas and water, whereas the lowest ETE values were computed for the 6CH 3 group both in gas and in water. These data clearly revealed the capability of TQ as a potential antioxidant molecule.

Frontier molecular orbitals (HOMO-LUMO)
Frontier molecular orbitals (FMO) or HOMO and LUMO denote the highest occupied molecular orbital and the lowest unoccupied molecular orbital respectively. FMOs are important parameters to predict the antioxidant activity of a compound [54]. The E HOMO represents the ability of electron donation, because it is energetically easiest to take away electron from this orbital, whereas E LUMO shows the capability to accept electron [55]. Molecules with lower E HOMO are more stable and less expected to give electrons, and the shape of HOMO determines the sites for free radical attack [48]. The calculated FMOs of TQ molecule in gas and water phases are represented in Fig. 3. FMO with traditional π-like molecular orbital appearance and distribution approximately over the entire molecule suggested antioxidant power of the molecule [41]. As shown in Fig. 3, HOMO and LUMO orbitals were distributed on an average over the entire TQ structure with little above the 5CH 3 and 6CH 3 groups, which directs their less involvement in antioxidant reactions. HOMO orbital was well-distributed over the 3CH and 12CH 3 groups, suggesting the probable contribution to antioxidant activity. The HOMO-1 orbital was mostly contributed to the π orbitals of double bonds and oxygen atoms. Theoretically, E HOMO is an effective indicator for free radical scavenging power [56]. LUMO orbital is mainly distributed over the ring with the two oxygen atoms. A small energy gap between LUMO and HOMO infers the softness, whereas a large gap indicates the hardness of the molecule. A soft molecule is more reactive than a hard one, and it can donate electrons more easily to the receiver. [57]. The global reactivity descriptors derived from E HOMO and E LUMO are calculated and presented in Table 5. The comfort of electron donation depends on the electronegativity. The values of hardness, softness, IP, EA, and electronegativity of TQ calculated in this study are very close to those computed for other molecules claimed as potential antioxidant [3].

Molecular electrostatic potential
Molecular electrostatic potential (MEP) is another supportive aspect for investigating and forecasting antioxidant capacity of a compound. The electron-deficient (positive) sites are the most likely sites for free radical attack [41]. It can be seen from the Fig. 4 that most electron-rich (negative) surfaces were positioned over the oxygen atoms, whereas the electrondeficient (positive) surfaces were located around the hydrogen atoms of the ring and the side chains in both gas phase and water phase. These upshots could be useful in the prediction of the positions for free radical attack of TQ molecule. However, the antioxidant reaction is an intricate and multifarious process and could be subjective to numerous aspects [58].

Thermodynamically favored mechanism
In order to investigate the antioxidant activity and mechanism of TQ, Gibbs free energy and enthalpy changes of the hydrogen abstraction (HA) from TQ to OH • radical were computed from the theoretically data of optimized geometry of the concerned reactants and products as per the simple principle of thermochemistry. The calculated ΔG and ΔH values for various CH groups of TQ against OH • radical are shown in Table 6. It is very interesting to report that all the ΔG and ΔH calculated herein were negative for the HA reaction which favor the antioxidant power of TQ through HAT mechanism.
The lowest ΔG and ΔH were observed for the 3CH group in all cases studied which is consistent with other thermochemical parameters such as BDE, PDE, and PA discussed in the earlier sections. In addition, from Table 6, it can easily be seen that

Conclusions
The antioxidant activity of TQ, a non-phenolic molecule, was successfully evaluated in gas and water phases using DFT calculation. The observed information revealed that the studied compound can pose as stronger antioxidant as other available ones. In both gas and water, 3CH and 12CH 3 positions were thermodynamically more feasible to scavenge free radical in comparison to 5CH 3 and 6CH 3 groups. In this regard, the HAT mechanism is preferable than SETPT and SPLET pathways. Obtained results clearly exhibit free radical scavenging potential of TQ molecule without possessing any phenolic group in its structure, which paves the way for other such compounds to be searched for their antioxidant power. This article has not been submitted elsewhere.  Code availability Not applicable.
Data availability All the data are provided in the manuscript and supplementary file.

Declarations
Consent to participate All the authors actively participated in this work.
Consent for publication If this manuscript is accepted for publication in this journal, we would not withdraw it.

Competing interests
The authors declare no competing interests.