Ligand-Protein Interactions : A Hybrid ab initio / Molecular Mechanics Computational Study

The enzymes Cyclooxygenase (COX) or prostaglandin-endoperoxide synthase (PTGS) are important in the synthesis of prostaglandins, which are the main mediating chemicals at inflammatory processes. The body produces two highly homologous COX isoforms, cyclooxygenase-1 (COX-1) and cyclooxygenase-2 (COX-2). COX-1 is involved in the production of prostaglandins which take part in physiological processes such as: protection of the gastric epithelium, maintenance of renal flow, platelet aggregation, neutrophil migration and also expressed in the vascular endothelium; Meanwhile COX-2 is inducible by proinflammatory stimuli. To counteract the symptoms of inflammation, nowadays is very frequent the use of nonsteroidal antiinflammatory drugs (NSAIDs); These drugs in addition to other benefits, can also cause side effects on people’s health (cardiovascular and respiratory problems, in the nervous system, among others). Due to the above, it is necessary to accelerate the investigations that allow to know in more detail the mechanisms of action that involve the use of natural plant products as pharmacological agents. In the present research, computational techniques are used, especially those based on the Molecular Docking Method to know the proteinligand interaction in systems of biological interest. It was possible to determine the structure-activity relationship in inflammatory processes involving the participation of a number of secondary plant metabolites such as luteolin, galangin, kaempferol, apigenin, morine and quercetin on the inactivation of the enzyme cyclooxygenase (COX-1 and COX -2). In the COX-1 / ligand coupling it was found that apigenin was the ligand which showed the lowest coupling energy with a value of -8.88, whereas quercetin showed the highest value (-6.65); for the coupling COX-2 / ligand the apigenin also showed the lowest energy value (-8.93), while kaempferol showed the highest value (-7.51). This shows that apigenin is the ligand with the best stability within the active site of both enzymes. On the other hand, quercetin showed the highest relative selectivity with a protein-ligand selectivity index (PLSI) of 1.17, while kaempferol showed an PLSI of 0.91. Taking into account the parameters of stability and selectivity we can say that of all the ligands used, quercetin would be the ideal to block COX-2.


Introduction
[5] In this work, we use a hybrid ab initio/mechanical-statistical technique to study the features of biological macromolecules, proteins, and their interactions with low molecular weight chemical species, called ligands.These protein-ligand systems, are a paradigmatic example of the structure-activity principle. 5,6We use here the Molecular Docking (MDock) technique, which allows us to predict the most appropriate conformation of a molecule (ligand) when it binds to another (protein) to form a stable complex. 68][9] Our goal is to understand the structure-activity relationship of the ligand-protein interactions from an atomistic point of view and, therefore, to attain relevant information to disclose the effect of a biological active molecule on a specific proteic target. 10e relevant proteic target, very useful in the clinical research field, is the Prostaglandin endoperoxide synthase (PTGS) or Cyclooxygenase (COX) enzyme. 11,12The first isoform is constitutive, and executes indispensable functions in the organism, while the second isoform intervenes in inflammatory processes, becoming the target of study for the search of novel anti-inflammatory compounds. 13th, COX-1 and COX-2, have a high structural similarity, only differing by the relative position of three amino acids.This structural similarity brings about non-selectivity between COX-1 and COX-2, e.g., non steroidal anti-inflammatory drugs (NSAIDs) are not selective for any form of COX, which causes serious consequences in the organism. 13,14e lack of NSAIDs-COX selectivity, has fueled the search for new bio-active molecules.Sometimes these species are obtained from natural sources, e.g., secondary metabolites from plants or animals, some of which are recognized to have a higher COX selectivity, helping to avoid or minimize the side effects produced by NSAIDs in the body. 15,16Within the secondary metabolites obtained from plants, we have a very wide group of chemicals known as flavonoids, which include: flavones, flavonols, antocyanidins and isoflavonols. 17These metabolites are recognized as anti-allergenic, antithrombotic, anti-oxidant and anti-inflammatory species. 18Some examples of vegetal species which are used to extract these compunds are: Dioclea grandiflora, Sophora flavescens, Swinglea glutinosa and Muntingia calabura. 15The aim of this work is to use tools of computational biochemistry to describe the structure-activity relationship between COX-ligand systems, which is of paramount importance in inflammatory inhibition, involving the action of secondary metabolites.

COX structural characterization
COX-1,2 enzymes belong to the oxidoreductase protein family, and forms complexes with selective inhibitors such as SC-558. 19We have used structural information from the PDB-database, as well as software such as PyMOL 20 to clean up the reported structures, for example, eliminating hemo groups, water molecules and ligands present in the structural data.

COX-1
From the Protein Data Bank (PDB) 21 , we have extracted 217 structures directly related with COX-1.
These structures were discriminated according to the species from which they have been obtained, and are distributed as follows: Homo sapiens (51), Bos taurus (29), Mus musculus (25), Ovis aries (22), Escherichia coli (12), Rattus norvegicus (12), Zea mays (10), others species (53). 21Currently, there is a record of 171 ligands that form complexes with COX-1, among which are flurbiprofen, arachidonic acid, ibuprofen, sodium ion and nitric acid. 21,22X-2 For COX-2, we have found 220 structures from the PDB, which were extracted from species like: Homo sapiens (51), Mus musculus (43), Bos taurus (28), Ovis aries (21), Escherichia coli (12), Rattus norvegicus (12), Saccharomyces cerevisiae (8), others (42). 21For this COX isoform, there are 178 ligands recognized, among which are the following: acetylsalicylic acid, flurbiprofen, urea, arachidonic acid and ibuprofen. 22r COX-1 we have used the model corresponding to Ovis aries (PDB ID: code 3N8Z) and for COX-2, the one for Mus musculus (PDB ID: code 1CX2).The experimental method used to dilucidate the tertiary structures of COXâ Ȃ Źs was reported to be Xray's crystallography. 19,23These structures were selected because, according to Kurumbail et al 19 , human COX-2 has a similarity of 87% with the murine enzyme, and the amino acid sequence that forms the active site is conserved between both species, therefore, we expect that human COX-2, and mainly the active site, will behave similarly to murine.

Molecular Docking
Before the MDock process, we add the hydrogen bonds missing to COX structures.This step is carried out because the hydrogen atoms have a small nuclei and are very mobile, so they are not detected by the X rays in crystallographic analysis.AutoDock Vina (AUV) 42 was used to guide the MDock Process.4][45] This software uses a semi-empiric force-field to evaluate the free-energy of conformations along the coupling simulations.At the same time, the force-field is parameterized using experimental information recover from numerous protein-inhibitor systems. 43We consider ten conformations for each protein-ligand pair.

Selectivity and Stability Indexes
We have developed an energetically based index of selectivity and stability for each COX-ligand pair.
The stability index is calculated by normalizing every Protein-Ligand interaction energy using the lowest calculated interaction energy.Likewise, the selectivity index was determined by dividing the binding energy of each ligand/COX-2 interaction between the binding energy of each ligand-COX-1 interaction, taking into account the correspondence between the ligands.Finally, the values obtained were normalized.

Molecular characterization of COX
The tertiary structures of COX are very similar, thus, they share the same folding units.Because of this, we use the COX-1 structure to characterize both enzymes.We can divide the tertiary structure of COX-1 in three folding units, showed in Figure 1.The first folding unit is formed by residues 34-72 in positions amino-terminal and form a small compact domain binded by three disulfide bonds, as Figure 2A shows.The global conformation of this domain, including the disposition of the disulfide bonds, is very similar to the conformation of the epidermal growth factor and binds covalently to the principal chain of the enzyme by another disulfide bridge, that is, Cys-37-Cys-159.The second unit is formed by residues 73-116 and form a structure of alpha helix A, B, C and D that are located at one side of the first folding unit, as we illustrate in Figure 2B.These alpha helix are highly amphipathic and represents a structural motif for the insertion of enzyme to the lipid bilayer.And the last folding unit consist of a large catalytic globular domain formed by residues 117-587.This catalytic domain is a globular structure that contains the cyclooxygenase and peroxidase active sites.Both sites, peroxidase and cyclooxygenase, comprise two differents lobes in the polypeptidic chain although interlaced with each other.The largest lobe is made by a structure formed by seven alpha helices, whereas the smallest lobe of the catalytic domain is formed by six alpha helices.These alpha helices form a set, more or less parallel, with its axes (see Figure 2C).Topological and molecular characterization of the active sites of COX

Characterization of the secondary structures of COX
The active sites of COX enzymes have a high similarity.Both active sites consist of a hydrophobic channel with an approximate lenght of 25 Angstroms.This channel originates in the Membrane Binding Domain (MBD) and extends to the nucleus of the globular domain. 19,23,46[49][50] The size of the active site of COX-2 is approximately about 20% more large that the COX-1 enzyme.This difference in the active sites of COX is very important due to the variation of three residues at positions 434, 513 and 523 of both enzymes, as we show in Figure 4.The variation in the size of both active sites has been an important feature for the development of specific NSAIDs for COX-2. 51Figure 5 shows, an superposition of the active sites of both enzymes with the residues that differentiate them, at the same point.

Geometry optimization
Table 3 shows the values obtained of minimization energies calculations and the evaluation of molecular HOMO and LUMO orbitals for each ligand.These values are ordered according to the HOMO energy from the best nucleophile to the best electrophile.The HOMO orbital is the last orbital that is doubly occupied and depicts the location of the pair of electrons most susceptible to electrophilic attack.Taking into account the values shown in Table 1, we can see that, of all the ligands, galangin is the best nucleophile, while apigenin is the worst nucleophile.On the other hand, Table 4 shows the geometry optimization energy values for each ligand in gaseous state, that is, before the coupling process occurs, as well as the optimization energies after the coupling process and the energy differences for COX-1/ligand interaction.Table 5 shows the values corresponding to COX-2/ligand interactions.

Molecular Docking
For each protein-ligand complex, we have selected the top ten conformations and for each one we evaluated the binding energy of the complex formed.Table 6 shows the interaction energy for each conformation for the Galangin/COX-1 system, as well as the amino acids and hydrogen bonds involved.For every conformation there is, at least, the formation of one hydrogen bond, except conformation 9 where it is shown the formation of 3 hydrogen bonds.The conformation number nine have the lowest energy-binding value, which means that in this conformation the ligand have more stability in the active site of enzyme compared with other conformations.
Table 6: COX-1/Galangin complex.1-10: obtained conformations and bindng energy values X: amino acids within the interaction orbit of the ligand, for each conformation.Line upwards: all aminoacids present in all of the conformations.Highlighted in green: amino acids interacting with the ligand and repeated in the indicated conformations.
for each COX-1/ligand pair and the amino acids that interact with the ligand in each conformation.
The interaction is remarkable by hydrogen bonds, attraction forces, electrostatic interactions and Van der Waals dispersion forces.Table 8 shows this information for COX-2/ligand system.Table 8 shows that conformation 7 (Apignenin) have the lowest value of coupling energy, namely, is the most stable conformation of the ligand in the active site of enzyme in the coupling process.
Galangin and Quercetin shows 2 conformations with the same amino acids and the same values of coupling energy.In this case, the ligands that follows Apigenin in high stability are: Morin, Quercetin, Luteolin, Galangin and Kaempferol, respectively.On the other hand, electrostatic terms has been the more relevants energetic terms to choose the conformations and the more stable ligands.Table 8 shows a comparison between the COX-1/ligand and COX-2/ligand interactions.Following the energetic terms, Apigenin have the most high stability in the active site of COX-1 with a value of -8.88 and a value of -8.93 in COX-2.While, the comparison of the coupling energy of Apigenin in COX-1 and COX-2, it shows that the stability is higher in the COX-1/Apigenin system.Stability Index SI We have been established a stability index according to the energy values (bind energies) of each protein-ligand interaction.We quantified the force exerted by the interaction of the ligand with the active site of enzyme.The interaction index is defined as follows: Where, E i represent the bind energy of each protein-ligand interaction and E minor is the minor bind energy of all protein-ligand interaction.The stability index is then normalized as follows: Where, NSI is the normalized stability index, E i represent the bind energy of each protein-ligand interaction, E minor is the minor bind energy of all protein-ligand interaction and E higher is the higher bind energy of all protein-ligand interactions.Quercetin have the lowest relative stability index of all ligands in complex with COX-1, while the same ligand has a higher stability index in complex with COX-2 (IE 0.18), which means that this ligand can be the structure with more energetic differences between both active sites.
On the other hand, Kaempferol have the higher stability index in COX-1 (IE 0.71) while this index is close to zero in COX-2.This suggest that Kaempferol is in a low energy state in the active site of COX-1.Table 10 shows a comparison of the energetic parameters and the formation of interactions between the ligands with IE values more representative.Figure 9 shows the interaction of Apigenin, Quercetin and Kamepferol with COX's, highlighting the relevant interacting AA's of active pockets and hydrogen bridges formed.

Ligand-Protein Selectivity Index of PLSI
The PLSI is a very important parameter to quantify the potential that some ligands have to inactivate a enzyme.This index shows the affinity for one or another isomeric form of COX in numerical and relative form.
The PLSI is calculated by dividing the binding energy of each COX-2/ligand interaction between the binding energy of each COX-1/ligand interaction, taking into account the correspondence between the ligands.
The PLSI is defined as follows: Quercetin was the ligand with a higher PLSI in the active sites of both enzymes, followed by Galangin, Morin, Apigenin, Luteolin and Kaempferol.
Figure 10 shows the PLSI.First, Quercetin have the higher value of selectivity with 1.17, while Kaempferol shows the lowest value with 0.91.Finally, we can notice that the difference in the PLSI values is not very marked among the ligands compared with the PLSI values, for which there are more marked differences.

Conclusions
We have developed a procedure to comparatively quantify the degree of selectivity and stability of chemical species in contact with a biologically active proteic agent.We have established a normalized stability index (NSI) based on docking binding protein-ligand energies, thus quantifying their degree of electrostatic (strong) and van der Waals (weak) intermolecular interaction.The NSI for every flavonoid in our set is defined as the normalized ratio between the calculated binding energy of each possible ligand-protein pair and the difference between the limit values.Also, from the ob-tained NSI's, we were able to define a protein-ligand stability factor (PLSF).The PLSF measures the affinity of a given set of bio-active chemical species towards a competing set of proteic targets.In our specific case, we have a set of COX-2/COX-1 inhibition competing flavonoids.Our chosen set of flavonoid species has been quantum-chemically characterized, highlighting its more elctrostatically active zones, as seen in the electronic density surface maps calculated from semiempirical hamiltonians.
Our results suggest that inside the active site of the COX's there is a complex interplay between the hydroxil groups of the ligand and the side chains of the sorrounding AA's.These hydroxil groups, are a critical variable, considering the set of flavonoids studied.However, the data shows that the formation of hydrogen bridges was not the most relevant factor for ligand stability/selectivity inside the active site.Thus, adscribing the potential stability or selectivity of a chemical species in terms of the number of hydrogen bridges formed with lateral chains of AA's, which are within the sphere of strong electrostatic or weak van der waals influence, is not the correct assesment.Although hydrogen bridges are important to characterize the set of acting ligand-protein forces, they do not determine the coupling stability.Based on the stability/selectivity criteria devised in this work, we could infer the independency between both factors.That means, a highly stable chemical species inside an protein active site does not guarantee an equally high degree of selectivity, wich could be unambiguously extracted from the apigenin data, which shows a high degreee of stability which is not equaled for the COX-2/COX-1 PLSI.From our flavonoid set of bioactive chemicals, we could mark Quercetin as the best candidate for the spcific selective inactivation of the highly homologous isomers of COX, having an intermediate level of stability and a large selectivity potential towards COX-2.
Finally, we are extending the present hybrid ab initio / mechanical statistics methodology to study wider sets of chemically bio-active relevant species related to activation/inactivation of protein targets, in order to make a more sensible screening of candidates to be subjected to in-vitro or in-vivo experimentation.

Figure 1 :
Figure 1: Folding units of COX-1.First unit is show in green, the second is show in orange and the third is show in violet.

Figure 2 :
Figure 2: Folding units of COX-1.(a) first folding unit in green and the three disulfide bonds are shown in yellow; (b) in orange is show the second folding unit and the four alpha helices (A, B, C and D); (c) in violet, the third folding unit is represent by alpha helices and beta sheets.

Figure 3 :
Figure 3: Amino acids of COX active site.(a) shown at the center of the figure are the amino acids that form the active site of the enzyme, highlighting the third folding unit.

Figure 4 :
Figure 4: Amino acids that differentiate the active sites of COX. A. Amino acid Ile434 is highlighting in yellow, His513 in green and Ile523 in red to COX-1.B. Amino acids of COX-2 are highlighting as: Val434 in yellow, Arg513 in green and Val523.
LUMO was carried out to know the distribution of the electrical charges of the molecules.Therefore, interactions with electrophilic groups are controlled by the HOMO orbital, while the nucleophilic attacks are controlled by the LUMO orbitals.Thus, Figure6allow us to observe the charge distribution and regions of the molecule with high electron density, that is, regions susceptible to electrophilic and nucleophilic attacks.These areas also show the electrostatic potential of the molecules, the red color shows the negative part and the blue color shows the positive part.The green and yellow regions show the poorest electron-rich zones.

Figure 6 :
Figure 6: Apigenin electronic density maps.(a) HOMO orbital: in red, the zone most susceptible to electrophilic attacks; (b) LUMO orbital: the region most susceptible to nucleophilic attacks.

Figure 7
Figure 7 shows, highlighted in red, five zones of high electronic density of Apigenin, in which the hydroxyl groups of the molecule are found.Due to the high electronic density, the formation of electrostatic interactions is possible.The amino acids that interact with Apigenin are: SER 530, ALA 378, ILE 377, PHE 381, PHE 209, PHE 205, VAL 228, LEU 534, GLY 533, ASN 375, PHE 529,

Figure 7 :
Figure 7: COX-2/Apigenin complex.The ligand is surrounded by amino acids that participate in the interaction with the active site of the enzyme.A hydrogen bond is shown in green.

Figure 8 (
Figure 8(a) represent the greater stability of Apigenin in the active site of COX-1, followed by Kaempferol, Luteolin, Morin, Galangin and Quercetin featuring the lowest stability.On the other hand, Figure 8(b) shows Apigenin as the more stable ligand in the active site of COX-2, followed by Morin, Quercetin, Luteolin, Galangin and Kaempferol, respectively.In the region of the ligands with the lowest stability, that is, 0.51, we find Luteolin, Morin, Galangin and Quercetin.

Table 3 :
Minimization energy values and molecular orbitals HOMO and LUMO [kJ/mol]It should be noted that the calculations of the energies of the molecular orbitals HOMO and

Table 5 :
Optimization and pos-coupling energy values of COX-2/ligand complex

Table 10 :
Energetic parameters and hydrogen bonds formation of Apigenin, Quercetin and Kaempferol.