Exploring interaction of TNF and orthopoxviral CrmB protein by surface plasmon resonance and free energy calculation.

Inhibition of the activity of the tumor necrosis factor (TNF) has become the main strategy for treating inflammatory diseases. The orthopoxvirus TNF-binding proteins can bind and efficiently neutralize TNF. To analyze the mechanisms of the interaction between human (hTNF) or mouse (mTNF) TNF and the cowpox virus N-terminal binding domain (TNFBD-CPXV), also the variola virus N-terminal binding domain (TNFBD-VARV) and to define the amino acids most importantly involved in the formation of complexes, computer models, derived from the X-ray structure of a homologous hTNF/TNFRII complex, were used together with experiments. The hTNF/TNFBD-CPXV, hTNF/TNFBD-VARV, mTNF/TNFBD-CPXV, and mTNF/TNFBD-VARV complexes were used in the molecular dynamics (MD) simulations and MM/GBSA free energy calculations. The complexes were ordered as hTNF/TNFBD-CPXV, hTNF/TNFBD-VARV, mTNF/TNFBD-CPXV and mTNF/TNFBD-VARV according to increase in the binding affinity. The calculations were in agreement with surface plasmon resonance (SPR) measurements of the binding constants. Key residues involved in complex formation were identified.

INTRODUCTION

Pro-inflammatory cytokines such as tumor necrosis factor (TNF), interleukins (IL-1, IL-3, IL-6, IL-8), interferons play the major role in the development of inflammation. TNF has received particular attention because its hyperproduction causes the development of chronic inflammations, including those of autoimmune nature. With this in mind, TNF is regarded as the principal target for novel biomedical technologies aimed at therapy of chronic pathologies such as rheumatoid arthritis, Crohn`s disease, among others. Different TNF inhibitors based on gene engineered products, which block the activity of this cytokine, are used in anti-cytokine therapy strategies. At present, a number of TNF-antagonists are available. Of these, Infliximab/Remicaid, Etanercept/Enbrall, and Adalimumab/Humira are applied for treating inflammatory diseases [1-4]. Infliximab is a chimeric monoclonal antibody against TNF containing the human constant and the mouse variable regions. It interacts with free and membrane-associated TNF with high affinity [1, 2]. Etanercept is a recombinant chimeric protein composed of two cystein-rich extracellular sub-domains of the human TNF receptor type II (hTNFRII) connected with an IgG1 Fc-fragment [3]. Adalimumab is a fully human monoclonal antibody against TNF [4].

The currently known inhibitors of TNF were efficient in the course of controlled studies of treatment of rheumatoid arthritis [5, 6]. Nevertheless, about 30-40% of the patients proved to be refractory (yielding poorly) in a real clinical setting; in less than a half complete or partial remission could be achieved whereas a third was compelled by secondary treatment inefficiency or adverse events 2-3 years after the treatment to give it up. The application of drugs acting as TNF inhibitors is restricted in a number of cases by tuberculosis [7], latent infections [8, 9]. The causes and underlying mechanisms of such complicating factors remain unclear. Obviously, there is a need for designing improved drugs of a new type that would be safer and block TNF activity more efficiently.

Another approach is now under consideration, the use of recombinant viral proteins that would block TNF activity. As known, certain viruses have acquired specific means during evolution for protecting themselves from the human immune system. In so doing, these viruses rely on the expression of proteins that can efficiently interfere with the organism’s protective response by blocking the activity of certain immune regulators. Examples are the orthopoxviruses that block TNF activity [10, 11]. Novel biomedical technologies for protecting humans rely on the use of these viral proteins. Such approaches are based on anti-cytokine action, and they are more advantageous in some characteristics compared to earlier developed drugs blocking the pathogenic role of TNF [12].

Modification of TNF-binding orthopoxviral proteins (TNFBPs) by site-directed mutagenesis is of importance because it increases TNFBPs binding efficiency and neutralizes TNF activity.Modeling of protein complexes structures and analysis of their interaction patterns would be then decisive. Computer molecular dynamics (MD) simulation is the most efficient approach for accomplishing this task.

Computer-based and theoretical methods have received growing recognition in analysis of the physical mechanisms of protein-protein complex formation. They extended the experimental data by providing new information about the nature of intermolecular interactions [13-25]. Recently, methods allowing estimation of free energy of the formation of molecular complexes are the first on the agenda. One is MM/GBSA, the Molecular Mechanics - Generalized Born Surface Area method [13-17, 26-29]. The MM/GBSA, together with MD simulation, gave valuable information about the structural and energetic interaction features by identifying the key amino acids involved in complex formation, exploring the contribution of particular mutations to the stability of a complex, thereby allowing delving deeper into the origin of chemotherapeutics [14-16]. Increased knowledge of the amino acids contributing to the stabilizing energy of the complexes (hot spots) would advance the development of potential drugs addressed to inhibition of the formation of protein-protein complexes through interaction with the hot spots.

We have performed a computer analysis of amino acid sequences of the viral TNF-binding CrmB proteins (TNF-BPs). These viruses included the variola virus (VARV), monkeypox virus (MPXV), and cowpox virus (CPXV), all pathogenic for humans, showing high homology and relatively small species-specific amino acid differences [30-32]. Individual recombinant proteins VARV-CrmB, MPXV-CrmB, and CPXV-CrmB have been synthesized in a baculovirus expression system in insect cells and isolated [33]. It has been shown that the CrmB proteins of VARV, MPXV, and CPXV differed in the efficiencies of inhibition of the cytotoxic effects of human (h), mouse (m) or rabbit TNFs for L929 mouse fibroblast cell line [33]. With this in mind, the aim of this work is computer MD simulation and comparative analysis of the mechanisms of molecular interaction of the hTNF or mTNF complexed with the N-terminal TNF-binding domains (TNFBD) of the CrmB proteins from CPXV and VARV.

The MD simulations of hTNF and mTNF homologous complexes with N-terminal TNFBDs of CrmB proteins allowed us to apply the MM/GBSA method in estimation of the formation free energy for these complexes. The results of the experimental analyses of the affinity for the formation of the complexes correlated well with the predicted estimates. Using additionally a MM-GBSA free energy decomposition protocol, we analyzed the contribution of individual amino acid residues to the stability of the complexes. It was found that the GLU138, TYR79, SER78, ALA25 amino acid residues in the hTNF and mTNF, also the LEU52, LEU56, LEU88, LEU89, ARG48 and ARG61 amino acids in the TNFBD-CPXV and TNFBD-VARV sequences are crucial in the stabilization of all the complexes. We predict that substitution of ASN63 by ASP63 results in a more efficient binding of the TNFBD-VARV to the human and mouse TNF compared with TNFBD-CPXV.

As a result, a better understanding of the molecular mechanisms of interaction of TNF-BPs of VARV or CPXV and TNF was gained. The results may be taken advantage of in designing site-specific mutagenesis experiments with the aim of increasing the affinity of interactions between TNF and the orthopoxviral TNF-BP and also of developing new potential therapeutic drugs.

MATERIAL AND METHODS

Homology Modeling

The 3D models of the protein structures were based on the X-ray structure of the TNFRII and mutant hTNF complex (pdbid:3ALQ). The 3D structures were generated using MODLLER 9v12 software. Sequence alignments were done using ClustalO (www.ebic.ac.uk/Tools/msa/clustalO). For each protein 50 models were generated. The models with lowest DOPE score were then subjected to energy minimization in vacuo using Amber12 with 15 steps of the conjugate gradient algorithm with distance-dependent dielectric constant ε = 4r. Verification of the 3D structure of the models was performed using PROCHEK [37, 38]. The constructed models were used for molecular dynamics simulations.

Molecular Dynamics (MD) Simulations and Binding Free Energy Calculations

The AMBER03 [39] force field as implemented in the AMBER12 [40] software package was used in the MD simulations. The complexes were solvated in a cube of TIP3P water with edges not smaller than 12 Å and cutoff of 12 Å for non-bonded interactions. The MD simulations were done with periodic boundary conditions. The particle mesh Ewald method [41, 42] was used for handling further electrostatic interactions. Time step was of 2 fs. The trajectory was written to file every 0.1 ps.

Each solvated complex was subjected to energy minimization using the conjugate algorithm [40]. Applying harmonic restraints with force constants of 2 kcal/(mol·A˚2) to all solute atoms, the canonic ensemble (NVT)-MD was generated for 50 ps during which the system was heated from 0 K to 300 K. The isothermal-isobaric (NPT)-MD was used during the next 50 ps of the MD simulation to adjust solvent density, following 1 ns equilibration without harmonic restraints. Next 15 ns (NPT)-MD simulation was carried out with gradually reducing distance restraints between CA atoms of CYS88 and SER116. It was done in order to stabilize structure of CRD3 B2 domain of TNF-BDs of CrmB. Afterwards 15 ns of unrestrained (NPT)-MD simulation was carried out before final 60 ns of production simulation (NPT)-MD. Bond lengths including those that involved hydrogen atoms were constrained using the SHAKE algorithm [43].

For the MM/GBSA method and per-residue free energy decomposition, snapshots were taken at 5 ps intervals from the corresponding final 60 ns of the MD trajectory. Explicit water molecules were removed from the snapshots. Energy was calculated using very large cutoff (999 Å). Binding free energy (ΔG) was calculated as

ΔGbind = ΔEgas + ΔGsolv - TΔSconf (1),

where ΔE is the interaction energy of molecules in the gas phase given as

ΔEgas = ΔEel + ΔEvdw (2),

where ΔE and ΔE represent electronic and van der Waals interactions, respectively. Solvation free energy (ΔG) was expressed as the sum of electrostatic solvation free energy (ΔG) and nonpolar contribution to solvation free energy (ΔG) was derived from the solvent accessible surface area (SASA) as ΔGΔ. S is the solvent-accessible molecular surface area calculated using a probe radius of 0.14 nm. ΔG and ΔG were calculated in a continuous solvent using the MM/GBSA model as implemented in the AMBER12 MM_PBSA.py software package. In calculation of the total free energy, the contribution of conformational entropy TΔS was disregarded because of high computational cost and low prediction accuracy. Moreover, we expect that it would not significantly affect relative free energy values, because of high similarity between investigated proteins [44]. Per-residue free energy decomposition was carried out using the protocol implemented in MM_PBSA.py [45].

Intermolecular contacts (hydrogen bonds and non-bonded contacts) were analyzed with DIMPLOT, which is part of the LIGPLOT software [46].

SPR-Analysis of Protein-Protein Interactions

A quantitative estimation of the efficiency of protein-protein interactions was obtained by using the ProteOn_XPR36 Protein Interaction Array System (Bio-Rad Laboratories, Inc.). The method is based on surface plasmon resonance. All experiments were carried out at 25 °С, flow rate was 30 µl/min. Recombinant proteins CrmB-VARV and CrmB-CPXV were covalently immobilized on chip GLM. GLM was activated for 2 min with a mixture of 1-(3-diethylaminopropyl)-3-ethylcarbodiimide hydrochloride (0.13 M) and N-hydroxysulfosuccinimide (0.03 M). Then, the viral proteins were immobilized at concentration of 50 µg/ml in 10 mM sodium acetate (pH 4.5) on chip surface, thereafter the surface was deactivated with 1M ethanolamine hydrochloride for 3 min, the captured immobilization level was 1000 RU. PBS phosphate buffered saline with 0.005% Tween 20 was used as running buffer. hTNF and mTNF (Invitrogen, USA) at 5 different concentrations ranging from 6.25 nM to 100 nM (at a flow rate of 30 µl/min) were used as analyte. After each measurement, the surface was regenerated with 10 mM glycine hydrochloride (pH 2.5). The binding interactions were fitted to a 1:1 Langmuir binding model. Measurements were done at the Orekhovich Institute of Biomedical Chemistry of the Russian Academy of Medical Sciences, Moscow.

Binding free energy of the two proteins was estimated using the formula ΔG = RT* ln (Kd), where R is the ideal gas constant of 1.985 cal/K∙mol, T is the absolute temperature 298 K, and Kd is the equilibrium affinity constant.

RESULTS AND DISCUSSION

The MODELLER 9v12 was applied to predict the 3D structures of the TNFBD-VARV, TNFBD-CPXV, hTNF and mTNF proteins [34-36]. To predict the 3D structure of the hTNF and mTNF, the mutant hTNF complex (PDB Id 3ALQ) was used as a template.

The amino acid sequences of hTNF and mTNF are 78% identical. The numbering of the amino acid sequences of mTNF and hTNF was based on the sequence of hTNF (Fig. ). The positions of hTNF SER85 (UNIPROT ID P01375) and mTNF SER88 (UNIPROT ID P06804) were taken as 1. The numbering of the amino acid sequences of TNFBD-VARV, TNFBD-CPXV was based on the sequence of TNFBD-CPXV (Fig. ).

To predict the 3D structure of the TNFBD-VARV and TNFBD-CPXV proteins, the structures of the human TNFBD of the hTNFRII (pdbId: 3ALQ) in complex with the mutant hTNF were taken as a template. Alignments of the amino acid sequences of TNFBD-VARV and TNFBD-CPXV with hTNFRII are given in (Fig. ). TNFBD of the hTNFRII showed over 40% identity to each of the examined viral TNFBDs, and it contained also a similar number of CRD domains and disulfide bonds participating in complex formation [47]. The TNFBD-CPXV and TNFBD-VARV sequences differ from each other by 6 amino acid residues and a deletion of 2 residues. For homology modeling we considered only TNF-binding domains of the CrmB-CPXV and CrmB-VARV proteins [33].

The hTNF/TNFBD-CPXV, hTNF/TNFBD-VARV, mTNF/TNFBD-CPXV, mTNF/TNFBD-VARV complexes were assembled by structural alignment of the models of the predicted 3D structures of hTNF, mTNF, TNFBD-CPXV, and TNFBD-VARV with the experimental 3D structure of the mutant hTNF/TNFBD-hTNFRII complex (PDB Id: 3ALQ).

More than 90% of all the amino acids were located at allowed regions of the Ramachandran plot and the common G-factor for each theoretical structure was close to the experimental structure (PDB Id:3ALQ). The predicted 3D structure of the mTNF/TNFBD-CPXV complex is shown in (Fig. , ).

Binding Free Energy

In the MM/GBSA calculations, stability of a system during simulation time is of great importance. To test the stability of the simulation, the root mean square deviation (RMSD) of Cα protein atoms over time was calculated. The RMSD was calculated relative to the starting complex structure in the production phase.

Figure shows that hTNF/TNFBD-CPXV, hTNF/ TNFBD-VARV, mTNF/TNFBD-CPXV and mTNF/ TNFBD-VARV simulated complexes were stable during the last 60 ns of MD simulation used for MM/GBSA free energy calculations. Deviation of the total energy, temperature and volume during the whole course of MD simulations was examined and shown to be stable (data not shown). As it can be seen from superimposed average structures (Fig. ), all complexes adapted similar conformations.

Table lists the contributions to binding free energy (ΔG) for the hTNF/TNFBD-CPXV, hTNF/TNFBD-VARV, mTNF/TNFBD-CPXV, mTNF/TNFBD-VARV complexes obtained by the MM/GBSA method. The components represent electrostatic and van der Waals interaction of proteins in the gas phase (ΔEΔEΔE), solvation free energy (ΔG) including polar (ΔG) and nonpolar (ΔG) contributions.

The calculated total binding free energies (ΔG) was in a good agreement with the experimental estimates (ΔGexp). Square of the Pearson`s correlation coefficient was equal R2=0.98 (Fig. ). It should be noted that the estimates calculated for the binding free energy did not reflect the absolute value of free energies of complex formation. This was because the contribution of conformational and translational entropy to the formation energy of the complexes was discounted, being computationally costly and giving very inaccurate predictions.

Since structures of the analyzed complexes showed a high similarity, it is expected that this assumption would not affect relative free energy values. Estimation of absolute free energies for which this term could be significant lies beyond this study.

As the data in Table show, binding affinity for the interaction between TNFBD-VARV and TNF is greater than that for TNFBD-CPXV, when hTNF and mTNF are considered separately. Concomitantly, it was found that mTNF could with greater affinity bind to the two viral TNF-binding proteins compared to hTNF.

A stabilizing effect of electrostatic interaction of the subunits was observed for all the complexes. An energetic desolvation penalty arose in solvent during complex formation. It resulted from loss of interaction of solvent and the contacting surface of the complex subunits. The balance between solvation energy and the energy of subunit interaction in the gas phase determined the stability of the formed complexes.

As can be seen from MM-GBSA free energy terms (Table ) the destabilizing contribution of polar solvent interactions exceeded the stabilizing contribution of the gas phase electrostatic interactions among the TNFBDs and TNF in the considered complexes. In turn, nonpolar Van-der-Waals interactions provided the main contribution to the stability of the complexes (Table ).

The higher affinity of interactions of TNF (human and mouse) with TNFBD-VARV in comparison to TNFBD-CPXV can be explained by different contribution of Van-der-Waals interactions in TNF/TNFBD-VARV and TNF/ TNFBD-CPXV complexes. One can see from Table , absolute gas phase electrostatic energy term grows up with the increase in total free energy of complexes formation, which indicates increase in a number of polar interactions with increase in binding affinity.

SPR-Analysis: Efficiency of the VARV-CrmB and CPXV-CrmB Recombinant Proteins Interaction with hTNF and mTNF

Production of recombinant virus proteins VARV-CrmB and CPXV-CrmB in Sf-21 insect cell culture and their purification were performed as described [48]. According to Alejo et al. [49], the CrmB protein is composed of two domains of which the N-terminal provides binding to TNF and the C-terminal the interaction with a set of chemokines. We have observed earlier that removal of the C-terminal chemokine-binding domain from the VARV-CrmB protein had no effect on its efficiency to inhibit TNF-induced cytotoxicity (in preparation).

The binding affinity of the recombinant viral receptors VARV-CrmB and CPXV-CrmB for the corresponding ligands was analyzed by SPR. The results showed that the Kd for the viral receptors differed depending on the species. Thus, K for the VARV-CrmB/hTNF complex was 2.48×10-9; K=3.62×10-10; K=4.10×10-9; K= 8.52×10-10. We have predicted earlier that the VARV-CrmB/mTNF is more stable than the VARV-CrmB/hTNF complex [44]. The SPR estimates of K for the VARV-CrmB/mTNF and VARV-CrmB/hTNF complexes were 1.14×10-4 and 7.59×10-4, respectively, thereby supporting the accuracy of the prediction. The SPR estimates of K for the CPXV-CrmB/mTNF and CPXV-CrmB/hTNF were 2.08×10-4 and 1.00×10-3, respectively.

Key Amino Acids Involved in the Formation of the Complexes

Per-residue free energy decomposition was done to specify the contribution of amino acid residues at the binding interface to the free energy of the complexes formation (Fig. ).

Up to 7 amino acid residues in hTNF and mTNF and up to 10 in the TNFBD-VARV and TNFBD-CPXV proteins making the greatest contribution to stabilization/destabilization of the complexes are distinguishable (|Δ G|>2 kcal/mol). Among them LEU52, LEU56, LEU88, LEU89, ARG48 and ARG61 in hTNF and mTNF, and ARG61, ARG48 in TNFBD-CPXV and TNFBD-VARV are key residues in stabilization of all the complexes according to the prediction. As can be seen, hydrophobic and polar amino acids have comparable contribution to total free energy. It should be noted that the above amino acid residues, the greatest contributors to energy formation of the complexes, are conserved in the TNF and TNF-binding protein families. This is further evidence of their functional significance. Figure shows that there are few polar amino acid residues which make a destabilizing contribution to interaction energy owing to high desolvation penalty (e.g. LYS57 in hTNF and mTNF, LYS97 in TNFBD-CPXV). However, their contribution is partly compensated by the contribution of the amino acids they interact with.

As was shown above TNFBD-CPXV and TNFBD-VARV have higher affinity to mouse TNF compared to human TNF. It can be explained by different contribution of amino acid residues to total free energy of complex formation in CRD1 region of TNFBD (Fig. ).

As seen (Table ), amino acid substitutions in the TNFBD-CPXV and TNFBD-VARV sequences do not contribute largely to the formation energy of the complex; they can, however, favor the creation of an environment appropriate for efficient interaction of the hot spot amino acids. Thus, ASN5/R->ASN5/R, ASN8/R->LYS8/R, ALA83/R->SER83/ R do not contribute directly to the stabilization of the complexes for the reason that they are not located at the region contacting with TNF, like the deletions of the ASN34/R and THR35/R amino acids are. The most difference between complexes of TNFBD-VARV and TNFBD-CPXV can be observed in the region around ARG61, LYS57 and LEU90 amino acid residues. It could be explained by amino acid substitution of ASP63->ASN63 leading to local conformational changes favorable for protein-protein interactions (Fig. ).

In the case of complexes with hTNF (Fig. , ) the substitution ASP63->ASN63 in TNFBD-VARV leads to formation of hydrogen bonds to phenolic OH-group of TYR107 (subunit A of hTNF) and loss of interaction with GLN59 (subunit A of hTNF). It also promotes interaction of GLN59 (subunit A of hTNF) with SER64 (TNFBD-VARV) leading to stabilization of the complex formation. In the case of complexes with mTNF (Fig. , ) substitution ASP63→ASN63 in TNFBD-VARV leads to complex structural rearrangements including perturbation of salt bridges.

CONCLUSION

Computer MD simulations together with laboratory experiments were used to analyze the mechanisms of interaction of the hTNF, mTNF and TNFBD-CPXV or TNFBD-VARV complexes, also to identify the most important amino acids involved in their formation. The structural models of the complexes were predicted using as a template the X-ray structure of the hTNF/TNFBD-TNFRII complex homologous with the four complexes. MD simulations were performed for models of the hTNF/TNFBD-CPXV, hTNF/ TNFBD-VARV, mTNF/TNFBD-CPXV, and mTNF/ TNFBD-VARV complexes. The MM/GBSA method was applied in free energy calculations. With this method, it was found that free energy of complex formation in which the mTNF was involved prevailed over that involving the hTNF. This was associated with a higher absolute free energy of the mTNF/TNFBD-VARV than that of the mTNF/TNFBD-CPXV complex, also with higher absolute free energy of the hTNF/TNFBD-VARV compared with the hTNF/TNFBD-CPXV.

Using the SPR, the efficiencies of the interaction of the recombinant TNF-BP-VARV and TNF-BP-CPXV with the hTNF or the mTNF were estimated. It was found that the experimental data for binding affinity were in agreement with predicted by MM/GBSA method.

The contributions of individual energy terms to free energy of complex formation calculated by the MM/GBSA method showed that polar interactions contributed significantly to the specificity of the interaction between TNF and the orthopoxvirus TNFBD. In all considered complexes, favorable free energy of complex formation was determined by nonpolar interactions.

The MM/GBSA free energy decomposition protocol enabled us to reveal the key amino acid residues involved in the stabilization of the protein complexes, also to explain the observed difference in their formation free energy. It was found that the amino acids GLU138, TYR79, SER78, ALA25 amino acids in the hTNF and the mTNF sequences, also the LEU52, LEU56, LEU88, LEU89, ARG48 and ARG61 in the TNFBD-CPXV and TNFBD-VARV sequen-ces are key residues in stabilization of all the complexes.

Based on analysis of MD trajectory, it may be concluded that substitution of ASN63 by ASP63 results in a more efficient binding of the TNFBD-VARV to the hTNF and the mTNF compared with TNFBD-CPXV.

The current work gave a deeper insight into the molecular mechanisms of the interaction of the TNF-binding protein of VARV and CPXV with TNF. In the future, we intend to plan site directed mutagenesis experiments to improve binding affinity of TNF and TNFBD using obtained computer models. The creation of anti-TNF agents based on the orthopoxviral TNF-binding protein is the ultimate goal.

CONFLICT OF INTEREST

The authors confirm that this article content has no conflict of interest.

Supplementary material is available on the publishers website along with the published article.

Table 1.

Contribution of Components to the Total Free Energy of Complex Formation.

	hTNF/TNFBD-CPXV	hTNF/TNFBD- VARV	mTNF/TNFBD-CPXV	mTNF/TNFBD- VARV
ΔEvdw	-105.2773	-109.58	-103.68	-122.18
ΔEel	-268.9	-668.7225	-1375.7	-1737.95
ΔEgas	-374.21	-778.3	-1479.38	-1860.13
ΔGGB	339.19	739.95	1431.76	1803.53
ΔGNP	-18.75	-19.89	-17.39	-16.97
ΔGsolv	320.44	720.06	1414.37	1786.56
ΔGel	70.29	71.22	56.06	65.58
ΔGbind	-53.7776	-58.24	-65.01	-73.56
ΔGexp	-11.58	-11.89	-12.53	-13.04

All values are given in kcal•mol-1. Components: ΔEvdv, van der Waals interaction energy; ΔEel, electrostatic interaction energy; ΔEgas, energy of complex formation in the gas phase; ΔGGB, polar interaction energy with solvent; ΔGNP, nonpolar interaction energy with solvent; ΔGsolv, solvation free energy; ΔGel, electrostatic interaction free energy; ΔGbind, total free energy of complex formation; ΔGexp, total free energy of complex formation, experimental.