Article pubs.acs.org/JPCB

Physical Nature of Fatty Acid Amide Hydrolase Interactions with Its Inhibitors: Testing a Simple Nonempirical Scoring Model Wiktoria Giedroyć-Piasecka,† Edyta Dyguda-Kazimierowicz,*,† Wiktor Beker,† Marco Mor,‡ Alessio Lodola,‡ and W. Andrzej Sokalski† †

Department of Chemistry, Wrocław University of Technology, Wrocław, Poland Pharmacy Department, Università di Parma, Parma, Italy

‡

S Supporting Information *

ABSTRACT: Fatty acid amide hydrolase (FAAH) is an enzyme responsible for the deactivating hydrolysis of fatty acid ethanolamide neuromodulators. FAAH inhibitors have gained considerable interest due to their possible application in the treatment of anxiety, inflammation, and pain. In the context of inhibitor design, the availability of reliable computational tools for predicting binding affinity is still a challenging task, and it is now well understood that empirical scoring functions have several limitations that in principle could be overcome by quantum mechanics. Herein, systematic ab initio analyses of FAAH interactions with a series of inhibitors belonging to the class of the N-alkylcarbamic acid aryl esters have been performed. In contrast to our earlier studies of other classes of enzyme−inhibitor complexes, reasonable correlation with experimental results required us to consider correlation effects along with electrostatic term. Therefore, the simplest comprehensive nonempirical model allowing for qualitative predictions of binding affinities for FAAH ligands consists of electrostatic multipole and second-order dispersion terms. Such a model has been validated against the relative stabilities of the benchmark S66 set of biomolecular complexes. As it does not involve parameters fitted to experimentally derived data, this model offers a unique opportunity for generally applicable inhibitor design and virtual screening.

■

thylketones, and acyl heterocycles.10 Other classes of inhibitors, characterized by an improved drug-like profile, have also been reported.11 These include piperazinyl-(pyridinyl)-urea and carbamate-based compounds,12 which have been shown to inhibit FAAH by covalently modifying the enzyme’s active site (i.e., through carbamoylation of the nucleophile Ser241 residue).13 Among these carbamoylating agents, N-alkylcarbamic acid aryl esters emerged as the first promising class of compounds capable of inhibiting FAAH in vivo, gaining considerable interest for the treatment of anxiety, inflammation, and pain.14 Clearly, the development of FAAH small-molecule inhibitors is a promising approach in the field of novel drug design, and many such attempts have been reported.12,15,16 Mor et al.12 have recently identified a lead compound with nanomolar activity, namely, N-cyclohexylcarbamic acid biphenyl-3-yl ester (compound URB524 in Figure 1). On the basis of this compound, a series of carbamic ester derivatives were designed and synthesized (scaffold on the right in Figure 1, variable group R1 indicates the substitution site), and the experimental activities of these compounds were measured.12 A number of quantum mechanics/molecular mechanics (QM/ MM) simulations allowed Mor and co-workers to establish that out of the two possible docking orientations, called positions I

INTRODUCTION Fatty acid amide hydrolase (FAAH, EC 3.5.1.99) is a membrane serine hydrolase responsible for the hydrolysis and inactivation of biologically active fatty acid ethanolamides,1−3 including endocannabinoid anandamide (arachidonoylethanolamine, AEA), and agonists of the peroxisome proliferatoractivated receptors, such as oleoylethanolamide and palmitoylethanolamide. FAAH is the only characterized mammalian member of the amidase signature (AS) family. This enzyme features a conserved serine- and glycine-rich stretch of approximately 130 amino acid residues, including a GXSXG motif, and a unique catalytic triad consisting of two serine residues and one lysine residue rather than the common serine−histidine−aspartate triad found in classical serine hydrolases.4 FAAH is also unusual among AS enzymes in that it is an integral membrane enzyme, whereas other known members of the AS family are soluble enzymes responsible for the hydrolysis of hydrophilic substrates.5−7 The pharmacological relevance of FAAH was realized when its chemical inactivation led to significant elevation in the endogenous concentration of AEA, potentiating the beneficial effects of AEA on pain and anxiety without eliciting the classical side effects of exogenous cannabinoid receptor agonists.2,8 FAAH therefore represents an attractive therapeutic target for the treatment of several central nervous system disorders.9 FAAH enzyme activity is blocked by a variety of classical serine hydrolase inhibitors, such as sulfonyl fluorides, fluorophosphonates, α-ketoesters, α-ketoamides, trifluorome© 2014 American Chemical Society

Received: June 14, 2014 Revised: November 21, 2014 Published: November 24, 2014 14727

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

be validated against benchmark data for the S66 set of biomolecular complexes.30

■

METHODS To investigate the physical basis of ligand binding within the FAAH active site, we have chosen seven amino acid residues (Met191, Ser193, Ser217, Ser376, Phe381, Phe432, and Val495) to represent the enzyme active site. The positions of these amino acid residues were obtained from FAAH complexes with docked ligands, which were provided by Mor et al. and described in refs 12 and 17. Figure 2 shows the FAAH active site representation with 4q (naphthalen-2-ylmethylcarbamic acid biphenyl-3-yl ester), which is one of the most potent FAAH ligands described in ref 12.

Figure 1. Chemical structure of lead compound URB524 and evidence of the key substituent at the nitrogen atom of the carbamate functional group.

and II, the latter promotes the inhibitory activity of the compounds of interest.17−19 This assumption was recently validated by X-ray crystallography,20 confirming that FAAH binds N-alkylcarbamic acid biphenyl-3-yl esters according to binding mode II. It is even more difficult to accurately determine binding energy to an enzyme target for covalent inhibitors because it depends on both stereoelectronic enzyme−inhibitor complementarity and the intrinsic chemical reactivity of the inhibitor.21 However, if the reactivity is kept essentially constant along a series of compounds, the inhibitory potency can be related linearly to the binding energy.22 In this scenario, it should be possible to predict differences in pIC50 for a series of inhibitors by simulating the enzyme−inhibitor recognition process (i.e., with molecular docking) and then estimating the binding affinity with a suitable and relatively accurate computational method. Such an approach can be undertaken due to the fact that even though the half-maximum inhibitory concentration (IC50) property formally reflects inhibitory activity (i.e., the effectiveness of a compound to inhibit enzyme activity), relative binding affinity can still be derived from IC50 values measured under consistent experimental conditions, which was the case herein.12 Following this reasoning, the binding energy differences calculated in our analysis were related to the experimental inhibitory activity. Statistical models correlating pIC50 values with classical (i.e., force field-derived) descriptors for FAAH inhibitors have been devised to provide a fair rationalization of experimental structure−activity relationship data;22 however, the first-principles understanding is still missing. The primary goal of this report is to contribute detailed nonempirical analyses on the physical nature of FAAH active site interactions with a reference class of inhibitors. The set of FAAH inhibitors employed herein includes members of the Nalkylcarbamic acid biphenyl-3-yl ester class with a variety of substituents at the nitrogen atom featuring different sizes, shapes, and electronic properties (R1 group indicated in Figure 1). Although currently available empirical scoring techniques used to estimate inhibitory activity in general provide an answer to the question of inhibitor affinity, there are systems in which these methods fail to reproduce experimentally observed trends.23−25 On the other hand, the ab initio theory of intermolecular interactions has been widely recognized as a well-established and accurate description of binding energy that could provide answers about not only the magnitude of the latter but also the physical nature of binding interactions.26−29 Systematic knowledge of the main factors contributing to receptor−ligand interactions could aid in the derivation of more computationally affordable yet sufficiently reliable methods of relating interaction energy to inhibitory activity. Therefore, in this study, we attempt to derive a simple nonempirical model based on electrostatic multipole and approximate dispersion terms. Additionally, such a model will

Figure 2. Structure of 4q inhibitor bound within the FAAH active site model.

The gaps that arose from cutting residues out of the protein structure were filled with hydrogen atoms optimized at the HF/ 6-31G(d) level of theory using Gaussian09.31 The set of inhibitors included in our analysis encompassed the URB524 structure (pIC50 = 7.20) together with 21 of its derivatives. R1 substituents characterizing the remaining inhibitors and the corresponding pIC50 values are listed in Table 1. As described in ref 12, the biphenyl scaffolds of all docked inhibitors were positioned almost identically (superimposition of the inhibitor poses is provided in the Supporting Information), suggesting that this particular fragment would contribute rather equally to the total binding energy. Therefore, it seemed reasonable to take into account only the R1 groups that differ between particular inhibitors. Accordingly, broken carbonyl carbon and ester oxygen bonds were filled with hydrogen atoms optimized at the HF/6-31G(d) level of theory. Such an approach has led to a reasonable estimate of the calculated interaction energy with a simultaneous decrease in computational time. The final model of FAAH−inhibitor complexes of the most potent 4q and weakest 4a inhibitor is given in the Supporting Information. Total FAAH−inhibitor binding energy was calculated in a pairwise fashion as the sum of the interaction energy components defined within the hybrid variational−perturbational theory (HVPT)32 for each amino acid residue from the FAAH active site model and a given inhibitor. Counterpoise correction was applied for the treatment of basis set superposition error.33 Interaction energy terms were calculated using the 6-31G(d) basis set34−36 with a modified37 version of the GAMESS38 program. Because the analyzed system is relatively large, application of the extensive basis set would compromise the efficiency of calculations. To test the reliability of the results obtained with the 6-31G(d) basis set applied herein, we compared the MP2/ 14728

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

obtained at the complete basis set (CBS) limit.30 Evaluation of the performance of the 6-31G(d) and the aug-cc-pVDZ39,40 basis sets in calculating the MP2 binding energy is given in the Supporting Information. Whereas the more extensive aug-ccpVDZ basis set provides a more accurate description of the MP2 interaction energy, as illustrated by the higher correlation coefficient value (R = 1.00 and 0.95 for MP2/aug-cc-pVDZ and MP2/6-31G(d), respectively), the correlation coefficient values for both basis sets indicate agreement between the results obtained with both 6-31G(d) and aug-cc-pVDZ basis sets and the reference CBS values. As shown in eq 1, the HVPT decomposition scheme applied herein allows us to partition the intermolecular interaction energy, EMP2, into the following contributions: the first-order electrostatic multipole term E(1) EL,MTP, representing the interaction of permanent multipole moments of two isolated molecules; the electrostatic penetration term E(1) EL,PEN, resulting from the overlap of charge distributions from interacting (1) monomers; the first-order exchange repulsion term EEX , originating from the Pauli exclusion principle; the delocalization (R) term EDEL , composed of the higher-order induction and exchange-induction terms; and the correlation term E(2) CORR. The latter could be further partitioned into the classical dispersion term EDas, resulting from intermolecular correlation, and intramolecular correlation contribution EINTRA−CORR to all electrostatic, exchange, delocalization, and dispersion terms.

Table 1. Structures of R1 Substituents and the Corresponding Experimental Affinitiesa of FAAH Inhibitors

(1) (1) (1) (R ) (2) EMP2 = E EL,MTP + E EL,PEN + E EX + E DEL + ECORR

(1)

In systems involving polar or charged molecules, the energy terms defined within the HVPT decomposition scheme are arranged according to the following hierarchy of gradually more complex and accurate levels of theory, which are obtained by incorporation of the subsequent corrections to the interaction energy.32 (1) (1) E EL,MTP < E EL < E(1) < ESCF < EMP2

(2)

As indicated below in eq 3, the consecutive approximations of interaction energy are related not only by strict hierarchy of

a

accuracy but also by increasing computational cost. In what follows, O(X) indicates scaling of the computational cost, where N and A stand for the number of atomic orbitals and atoms, respectively. The simplest term (i.e., the electrostatic multipole (1) component EEL,MTP ) is calculated using atomic multipole expansion. The penetration term is obtained as the difference between the first-order perturbation expression E(1) EL and the

Inhibitory activity values taken from ref 12.

6-31G(d) interaction energy values evaluated for the benchmark S66 set of biomolecular complexes to the reference results 14729

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

multipole contribution E(1) EL,MTP. The exchange repulsion term, , is obtained as the difference between the first-order E(1) EX Heitler−London energy, E(1), and the electrostatic term, E(1) EL . Higher-order delocalization energy is defined as the difference between the counterpoise-corrected self-consistent field variational energy, ESCF, and E(1) energy. Finally, the correlation contribution, E(2) CORR, is obtained as the difference between the second-order Möller−Plesset interaction energy, EMP2, and ESCF energy. As mentioned, the introduction of consecutive corrections results in not only improved precision of the corresponding energy values but also increased computational cost. HVPT is a valid method to achieve a reasonable compromise between the expected accuracy of calculations and the computational effort required by the latter, rendering this approach applicable even for relatively large systems, such as enzyme−inhibitor,32,41,42 enzyme−reactant,43,44 and RNA−intercalator complexes.45,46 In most of these systems, predictions based on the electrostatic term correlated with the most complete MP2 model. Moreover, relative stabilities based on the electrostatic term appeared to be less sensitive to inaccurate structural data47,48 due to elimination of the repulsive exchange term. Considering the prevalence of neutral, or even nonpolar, amino acid residues in the FAAH active site,49 accompanied by neutral FAAH inhibitors, the electrostatic energy term might be less indicative in terms of predicting inhibitory activity because the contribution of computationally expensive dispersive forces might be non-negligible. The dispersion corrections have mostly been developed in conjunction with density functional theory methods which, in principle, neglect dispersive van der Waals interactions.50−53 A promising approach to calculate the dispersion approximation, EDas, has recently been proposed54 and appears to be well-suited for describing noncovalent interactions. In particular, it was shown to perform exceptionally well describing the interaction energy within a set of hydrogen-bonded alcohol dimers in which both electrostatic and dispersion terms appear to play an important role.55 EDas contribution can be estimated at a very low computational cost from a universal function approximating dispersion and exchange-dispersion energies for a training set of dimers.54 The performance of various scoring methods could also be compared in terms of the success rate (Npred) of predicting relative affinities. Using nonparametric statistics, Npred is defined as the percentage of concordant pairs of complexes among all pairs of complexes within a given set.47,48 As a reference, we used the most accurate ΔECBS CCSD(T) interaction energy values extrapolated to the complete basis set at equilibrium geometry for each dimer from the S66 set.30 It is remarkable that the costly MP2 results fail to reproduce relative stabilities for distances shorter than equilibrium, whereas much simpler models composed of electrostatic multipole term exhibit reasonable agreement within this range of intermolecular separation.48 This could be extremely significant in a typical situation in which structures of enzyme−ligand complexes are generated by approximate docking procedures, which frequently result in shortened artifactual contacts.42 The results presented in Figure 3 indicate that, for a range of monomer separations, the most cost-effective electrostatic multipole model, supplemented by an approximate dispersion term, indeed exhibits improved predictivity at both the near equilibrium geometry and larger distances.

Figure 3. Success rate (Npred) of predicting the relative stability of 2145 dimer pairs within the S66 set of 66 typical biomolecular complexes30 as a function of the R/Req factor scaling the shortest (1) + EDas denotes an electrostatic intermolecular contact. EEL,MTP (2) multipole and approximate dispersion term, and E(1) EL,MTP + ECORR stands for an electrostatic multipole term augmented by complete correlation contribution.

■

RESULTS AND DISCUSSION A list of the inhibitor compounds considered herein, together with their FAAH inhibitory potency, is reported in Table 1. Table 2 presents the values of electrostatic multipole and Table 2. Contributions to FAAH−Inhibitor Interaction Energya

a

inhibitor

E(1) EL,MTP

E(1) EL,PEN

E(1) EX

E(R) DEL

E(2) CORR

EDas

4q 4u 4v 4r 4s 4h 4t 4f 4g 4e 4p 2 4i 4d 4l 4m 4o 4n 4c 4b 4k 4a

−5.1 −7.0 −5.3 −5.4 −7.2 −5.1 −6.6 −5.8 −5.2 −5.8 −6.0 −5.6 −6.4 −6.6 −5.9 −4.1 −3.7 −1.5 −5.0 −5.8 −2.6 −5.3

−4.0 −4.4 −5.2 −4.3 −6.8 −3.8 −7.0 −4.9 −2.0 −3.8 −5.0 −3.9 −6.0 −4.3 −3.9 −3.3 −4.1 −2.3 −2.4 −4.1 −4.6 −4.5

12.7 15.1 16.2 14.6 19.9 11.5 21.8 15.3 6.1 11.8 16.1 12.8 19.5 15.0 12.9 10.4 12.6 6.5 8.4 13.0 12.2 14.0

−3.3 −3.7 −3.6 −3.5 −4.4 −2.9 −4.9 −3.8 −1.9 −3.4 −4.3 −3.0 −5.5 −4.2 −3.5 −2.9 −3.2 −1.4 −1.9 −3.9 −1.8 −3.8

−5.2 −6.6 −8.5 −6.1 −7.7 −3.3 −7.3 −4.7 −2.6 −3.2 −5.5 −4.1 −3.6 −3.1 −3.6 −2.9 −5.4 −4.4 −2.9 −1.8 −4.6 −1.7

−11.6 −14.4 −17.2 −13.6 −16.4 −9.6 −17.7 −12.6 −7.6 −9.8 −13.2 −10.9 −11.0 −10.0 −9.8 −8.4 −12.4 −9.4 −8.3 −8.0 −10.9 −7.8

In units of kcal mol−1.

penetration, exchange, delocalization, correlation, and dispersion contributions to the total interaction energy. The approximate dispersion term EDas has been introduced as a computationally affordable alternative to the complete (2) correlation energy E(2) CORR. Whereas the absolute ECORR and EDas values exhibit pronounced differences, these two terms coincide with each other to a great extent as their mutual 14730

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

Table 3. Total FAAH−Inhibitor Interaction Energya at Consecutive Levels of Theory inhibitor

pIC50b

E(1) EL,MTP

E(1) EL

E(1)

ESCF

EMP2

E(1) EL,MTP + EDas

4q 4u 4v 4r 4s 4h 4t 4f 4g 4e 4p 2 4i 4d 4l 4m 4o 4n 4c 4b 4k 4a

8.28 8.27 8.11 8.03 7.89 7.47 7.40 7.28 7.27 7.24 7.23 7.20 7.18 6.95 6.86 6.76 6.67 6.31 6.28 6.16 5.39 4.86

−5.1 −7.0 −5.3 −5.4 −7.2 −5.1 −6.6 −5.8 −5.2 −5.8 −6.0 −5.6 −6.4 −6.6 −5.9 −4.1 −3.7 −1.5 −5.0 −5.8 −2.6 −5.3 −0.49 0.3

−9.1 −11.4 −10.5 −9.7 −14.0 −8.9 −13.6 −10.7 −7.2 −9.6 −11.0 −9.5 −12.4 −10.9 −9.8 −7.4 −7.8 −3.8 −7.4 −9.9 −7.2 −9.8 −0.43 0.5

3.6 3.7 5.7 4.9 5.9 2.6 8.2 4.6 −1.1 2.2 5.1 3.3 7.1 4.1 3.1 3.0 4.8 2.7 1.0 3.1 5.0 4.2 0.15 0.5

0.3 0.0 2.1 1.4 1.5 −0.3 3.3 0.8 −3.0 −1.2 0.8 0.3 1.6 −0.1 −0.4 0.1 1.6 1.3 −0.9 −0.8 3.2 0.4 −0.01 0.3

−4.9 −6.6 −6.4 −4.7 −6.2 −3.6 −4.0 −3.9 −5.6 −4.4 −4.7 −3.8 −2.0 −3.2 −4.0 −2.8 −3.8 −3.1 −3.8 −2.6 −1.4 −1.3 −0.83 0.3

−16.7 −21.4 −22.5 −19.0 −23.6 −14.7 −24.3 −18.4 −12.8 −15.6 −19.2 −16.5 −17.4 −16.6 −15.7 −12.5 −16.1 −10.9 −13.3 −13.8 −13.5 −13.1 −0.67 0.9

Rc SEd

In units of kcal mol−1. bInhibitory activity values taken from ref 12. cCorrelation coefficient between the energy obtained at a given level of theory and the experimental inhibitory activity. dStandard error (kcal mol−1).

a

correlation coefficient is equal to 0.96 (see corresponding plot in the Supporting Information). Apparently, intermolecular correlation is the main factor contributing to changes in E(2) CORR. The total interaction energy for all inhibitors calculated at each level of theory is given in Table 3. As demonstrated by a significant value of the correlation coefficient (R, Table 3), MP2 binding energy correlates with the experimental inhibitory activity. Compared to the MP2 correlation coefficient, R values for the more approximate levels of theory are substantially worse. This indicates an essential role of correlation effects, especially the dispersion term EDas, presumably arising from the absence of polar or charged residues in the variable part of the inhibitor structures. Noticeably, the remaining interaction energy components cancel out to a significant degree. As a result, only the most robust MP2 level of theory provides a satisfactory value of the correlation coefficient R, suggesting that, due to the neutral nature of the FAAH active site, the proper model of FAAH inhibitory activity requires that all intermolecular binding-force contributions are included. The composition of the MP2 interaction energy for each inhibitor is shown in Figure 4, which provides the binding contribution of each amino acid residue. Clearly, the methionine residue Met191 has a major influence on total binding affinity. Nevertheless, exclusion of any other residue results in decreased values of the correlation coefficient, indicating that all seven residues are essential to properly represent the specificity of FAAH inhibitor binding. Interestingly, within the residues constituting the FAAH catalytic triad as reported in the literature (Ser241−Ser217−Lys142),6,7,12 only Ser217 contributes significantly to the experimental binding affinity.

Figure 4. Contributions of selected amino acid residues to FAAH inhibitor binding.

The following equation calculates the predicted inhibitory activity from the computed MP2 binding energy values in kcal mol−1 pIC50 =

5.813 − EMP2 1.384

(4)

Overall, reasonable agreement was obtained between values of the negative logarithm of the theoretical binding affinities as derived from eq 4 and their experimental equivalents (Figure 5). These results indicate that once the correlation effects are 14731

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

docking and/or optimization procedures,42 the E(1) EL,MTP term appears to provide a more suitable description of the binding characteristics.47,48 Thus, application of the electrostatic multipole-based interaction energy estimates could be beneficial on the basis of their relative insensitivity to inaccuracy in structural data and elimination of the excessively large repulsion term. Analysis of the data presented in Table 4 suggests that for inhibitor-bound complexes of phenylalanine ammonia-lyase (PAL),41 leucine aminopeptidase (LAP),41 or urokinase-like plasminogen activator (uPA)42 (i.e., systems encompassing charged residues or receptor−ligand hydrogen bonds), satisfactory correlation with the experimental data could be (1) or obtained at the most approximate electrostatic EEL levels of theory. On the contrary, electrostatic multipole E(1) EL,MTP both of these models failed to yield correlation for FAAH inhibitors. Such an unusual situation might arise from the fact that the FAAH active site is composed mostly of nonpolar amino acid residues interacting with neutral inhibitors. Accordingly, the overall effect on the electrostatic energy is moderate, indicating that the dispersive forces, which originate from interactions with polarizable side chains, might be the key contribution to the FAAH inhibitory effect of the molecules (1) under study. Indeed, adding the E(2) CORR term to EEL,MTP results in an R value of −0.74, which constitutes considerable improvement in terms of agreement with the experimental data. Unfortunately, calculation of the correlation energy E(2) CORR or its main component (i.e., dispersion interaction) is extremely costly because it scales as at least O(N5) and, additionally, requires the use of very extended basis sets with polarization functions. Pernal et al. have derived a general formula approximating the accurate dispersion interaction energy, EDas, in the form of atom−atom potentials, including damping terms.54 Recently, this formula has been successfully tested to reproduce the binding energy values of alcohols.55 In addition, it improves the success rate of predicting relative stabilities in the full S66 test set of typical biomolecular complexes (Figure 3). In terms of the correlation with experimental binding affinity (2) values, the E(1) EL,MTP+ ECORR model yields substantially improved predictions (R = −0.74) over those obtained for E(1) EL,MTP alone (R = −0.49). However, the E(1) EL,MTP + EDas model described by R = −0.67 and scaling with O(A2) constitutes the best (2) compromise, as both the EMP2 and E(1) EL,MTP + ECORR models are

accounted for (the dispersion term in particular), the potency of a given FAAH inhibitor might be predicted.

Figure 5. Inhibitory activities predicted from MP2 interaction energy as a function of experimental binding affinity.

Figure 5 demonstrates that inhibitors 2 and 4i, which have essentially the same inhibitory potency (Table 1), differ in terms of their MP2 binding energy. Despite these two compounds having similar R1 substituents, the cyclohexanyl group does not occupy the same spatial position, probably because of the presence of an additional methylene group in inhibitor 4i (the corresponding spatial arrangement is presented in the Supporting Information). As indicated by comparing the shortest distances between a given inhibitor and the FAAH amino acid residues (Supporting Information), the Met191 residue is closer to the 4i inhibitor (the minimum distances are equal to 2.0 and 1.7 Å for compounds 2 and 4i, respectively). In the case of inhibitor 2, the overall MP2 binding energy is mostly the result of Met191 contribution, which amounts to −3.1 kcal mol−1 (Figure 4 and Table S1 in the Supporting Information). The stabilizing influence of the latter does not occur in the FAAH−4i inhibitor complex, probably due to excessive shortening of the Met191−inhibitor distance, as suggested by a significant increase in the exchange repulsion term and the resulting destabilization of the Met191− inhibitor complex observed at the E(1) level of theory (Supporting Information). In the case of suboptimal intermonomer distances that might be obtained with approximate

Table 4. Comparison of Correlation Coefficient R and Standard Error Values Associated with Nonempirical Models Obtained from HVPT Analysis Performed for Inhibitors of Several Enzymes method

PALa

LAPb

EMP2 ESCF E(1) HL E(1) EL E(1) EL,MTP number of hydrogen bond donors and acceptors number of charged residues number of inhibitors in the correlation

−0.99 (5.0) −0.99 (5.1) −0.92 (4.6) −0.88 (6.4) −0.72 (10.7) 4 1 6 d

−0.97 −0.95 −0.94 −0.93 −0.96

(24.5) (28.7) (101.8) (92.2) (30.3) 4 2f 5

uPAc

FAAH

−0.79 (3.2) −0.58 (2.5) 0.87 (6.2) −0.94 (15.8) −0.99 (10.2) 3 1 5

−0.83 (0.3) −0.01 (0.3) 0.15 (0.5) −0.43 (0.5) −0.49 (0.3) 0e 0 22

a

Phenylalanine ammonia-lyase (PAL) results taken from ref 41. bLeucine aminopeptidase (LAP) results taken from ref 41. cUrokinase-like plasminogen activator (uPA) results taken from ref 42. dStandard error values in parentheses (kcal mol−1). eRefers to the fragments of FAAH inhibitors accounted for in the calculations. fIn addition to two charged amino acid residues, two zinc cations (Zn2+) are also present in the LAP active site. 14732

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

much more costly (scaling with O(N5), see eq 3) and therefore impractical. It must be kept in mind that the enzyme−ligand (or inhibitor, in particular) interaction is an extremely complex phenomenon.56 In particular, the process of ligand binding involves many other factors besides nonbonded intermolecular interactions, for instance, desolvation effects (removal of water molecules from an enzyme active site upon ligand binding) and residual flexibility (entropic stabilization of the protein).57−59 Nonempirical calculations accounting for the effects mentioned above would be extremely costly. However, assuming that the most important component of the protein−ligand interaction is the electronic part of the enthalpic contribution on the basis of its prevalence, and that any other enthalpic terms as well as the entire entropic term can be neglected in the binding affinity calculations,60 allows for a good approximation of ligand affinity within a reasonable amount of time and without introducing any empirical parameters. According to a recent study by Rao et al., solvation and entropic effects can be neglected without introducing significant errors in the calculations.61 As a result, an ab initio approach, dependent only on a set of fundamental constants and the results of high-level quantum chemical calculations, becomes a universal method for studying molecular interactions. This is an advantage over any method employing empirical parameters derived from an arbitrarily chosen set of experimental data, wherein the results can be valid for one system but not necessarily another; success depends on the similarity of the system under investigation to the training set of molecules that were used to derive the parameters implemented within a given method.62−65 Thus, ab initio methods should ideally be incorporated into the process of ligand design. The magnitude of the standard error values accompanying the respective enzyme−inhibitor systems compared in Table 4 warrants additional explanation. Depending on the prevalent type of interactions and the size of the system under study, the calculated binding energy values might differ significantly between particular inhibitors (e.g., in the case of the LAP enzyme featuring a binding site composed of two zinc ions and two other charged residues), or they might exhibit less pronounced differences (FAAH results, discussed herein). Another factor affecting the standard error value is the size of the sample (i.e., the number of inhibitors included in a given set; see Table 4 for comparison). Because these determinants are different in the cases summarized in Table 4, a direct comparison of the standard error values characterizing particular enzyme−inhibitor systems is not valid. To compare nonempirical results with available empirical scoring functions, we used Discovery Studio 3.5 Suite66 to calculate binding, electrostatic, and entropic energies.67 Scoring methods available in Discovery Studio 3.5 (i.e., LigScore,68 Picewise Linear Potential (PLP),69,70 Jain,71 Potential of Mean Force (PMF),72,73 and Ludi74,75 functions) were employed in the analysis of FAAH inhibitor binding. The performance of ab initio methods is further validated by analysis of correlation coefficient values and the percentage of successful relative-stability predictions obtained with empirical scoring methods (Table 5). Except for PMF and Ludi, the remaining empirical functions do not correlate with the experimental binding affinity. The greater discriminatory capabilities of the ab initio model over the predictions of inhibitory potency obtained with empirical scoring functions indicate the possibility of using the more robust QM

Table 5. Performance of ab Initio and Empirical Scoring Methods for Reproducing the Activity of FAAH Inhibitors method

Ra

Npred (%)b

EMP2 (2) E(1) EL,MTP + ECORR E(1) + E EL,MTP Das PMF Ludi1 PLP2 Jain LigScore1

−0.83 −0.74 −0.67 −0.72 −0.62 −0.51 −0.48 +0.25

83.1 78.4 74.9 77.1 73.2 65.8 71.4 44.6

a

Correlation coefficient between the energy obtained with a given method and the experimental inhibitory activity. In the case of empirical scoring functions, for which a higher score indicates greater inhibitory potency, the opposite of the correlation coefficient value is given to allow for direct comparison with the results of nonempirical binding energy calculations, wherein the better inhibitor is associated with a lower interaction energy value. bPercentage of successful predictions.

calculations throughout the virtual screening process, especially regarding the computationally affordable E(1) EL,MTP + EDas model. The computational cost of the latter is as favorable as that of commonly applied scoring functions (see the Supporting Information for a detailed comparison). The advantage of using a nonempirical approach is the lack of calibration and training on experimentally determined affinity data. Undoubtedly, empirical scoring functions fitted to reproduce some selected experimental results may work very well for some types of ligands and protein families, but they lack predictive abilities in many other instances.65 The general utility of such a simple nonempirical model will be tested in our lab against more protein−ligand complexes.

■

CONCLUSIONS Docking has become a standard technology in drug discovery to virtually screen libraries of existing or hypothetical compounds with the goal of identifying new active ligands that can bind to a target protein and then to predict the binding modes and affinities of those ligands.76,77 Although docking methods are quite efficient tools for these studies, the reliability of the docking predictions is currently limited by several methodological shortcomings, including the energetic quantification of protein−ligand binding.63,64 In this study, we used ab initio calculations to derive a scoring approach to rank a series of docked FAAH inhibitors featuring highly lipophilic substituents. Relative ab initio binding energy values for a series of FAAH inhibitors correlate with their respective experimental affinities, indicating that the activity of a potential FAAH inhibitor could be predicted within a reasonable error at the MP2 level of theory. Because of the mostly nonpolar character of the variable R1 substituents and the considerable influence of correlation effects, none of the simpler levels of theory seem to be sufficient for obtaining satisfactory results, although the proposed E(1) EL,MTP + EDas model could serve as an inexpensive alternative to EMP2 calculations. The counterintuitive importance of electrostatic interactions could best be illustrated by the dispersion-dominated stacking interactions in nucleic acid bases, wherein the much smaller electrostatic term is the best predictor of the relative stabilities.78 On the other hand, our current results indicate 14733

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

degradation modulates the expression of inflammatory mediators in the hypothalamus following an immunological stressor. Neuroscience 2012, 204, 53−63. (4) McKinney, M.; Cravatt, B. Structure and function of fatty acid amide hydrolase. Annu. Rev. Biochem. 2005, 74, 411−432. (5) Bracey, M.; Hanson, M.; Masuda, K.; Stevens, R.; Cravatt, B. Structural adaptations in a membrane enzyme that terminates endocannabinoid signaling. Science 2002, 298, 1793−1796. (6) Labar, G.; Michaux, C. Fatty acid amide hydrolase: from characterization to therapeutics. Chem. Biodiversity 2007, 14, 1882− 1902. (7) Otrubova, K.; Ezzili, C.; Boger, D. The discovery and development of inhibitors of fatty acid amide hydrolase (FAAH). Bioorg. Med. Chem. Lett. 2011, 21, 4674−4685. (8) Hauer, D.; Weis, F.; Campolongo, P.; Schopp, M.; BeirasFernandez, A.; Strewe, C.; Giehl, M.; Toth, R.; Kilger, E.; Schelling, G. Glucocorticoid-endocannabinoid interaction in cardiac surgical patients: relationship to early cognitive dysfunction and late depression. Rev. Neurosci. (Berlin, Ger.) 2012, 23, 681−690. (9) Petrosino, S.; Marzo, V. D. FAAH and MAGL inhibitors: therapeutic opportunities from regulating endocannabinoid levels. Curr. Opin. Invest. Drugs (BioMed Cent.) 2010, 11, 51−62. (10) Seierstad, M.; Breitenbucher, J. Discovery and development of fatty acid amide hydrolase (FAAH) inhibitors. J. Med. Chem. 2008, 51, 7327−7343. (11) Minkkila, A.; Saario, S.; Nevalainen, T. Discovery and development of endocannabinoid-hydrolyzing enzyme inhibitors. Curr. Top. Med. Chem. (Sharjah, United Arab Emirates) 2010, 10, 828−858. (12) Mor, M.; Lodola, A.; Rivara, S.; Vacondio, F.; Duranti, A.; Tontini, A.; Sanchini, S.; Piersanti, G.; Clapper, J.; King, A.; et al. Synthesis and quantitative structure−activity relationship of fatty acid amide hydrolase inhibitors: modulation at the N-portion of biphenyl3-yl alkylcarbamates. J. Med. Chem. 2008, 51, 3487−3498. (13) Alexander, J.; Cravatt, B. Mechanism of carbamate inactivation of FAAH: implications for the design of covalent inhibitors and in vivo functional probes for enzymes. Chem. Biol. 2005, 12, 1179−1187. (14) Kathuria, S.; Gaetani, S.; Fegley, D.; Valino, F.; Duranti, A.; Tontini, A.; Mor, M.; Tarzia, G.; Rana, G. L.; Calignano, A.; et al. Modulation of anxiety through blockade of anandamide hydrolysis. Nat. Med. (N.Y., NY, U.S.) 2003, 9, 76−81. (15) Ahn, K.; Johnson, D.; Mileni, M.; Beidler, D.; Long, J.; McKinney, M.; Weerapana, E.; Sadagopan, N.; Liimatta, M.; Smith, S.; et al. Discovery and characterization of a highly selective FAAH inhibitor that reduces inflammatory pain. Chem. Biol. 2009, 16, 411− 420. (16) Butini, S.; Brindisi, M.; Gemma, S.; Minetti, P.; Cabri, W.; Gallo, G.; Vincenti, S.; Talamonti, E.; Borsini, F.; Caprioli, A.; et al. Discovery of potent inhibitors of human and mouse fatty acid amide hydrolase. J. Med. Chem. 2012, 55, 6898−6915. (17) Lodola, A.; Mor, M.; Rivara, S.; Christov, C.; Tarzia, G.; Piomelli, D.; Mulholland, A. Identification of productive inhibitor binding orientation in fatty acid amide hydrolase (FAAH) by QM/ MM mechanistic modelling. Chem. Commun. 2008, 2, 214−216. (18) Lodola, A.; Capoferri, L.; Rivara, S.; Chudyk, E.; Sirirak, J.; Dyguda- Kazimierowicz, E.; Sokalski, W.; Mileni, M.; Tarzia, G.; Piomelli, D.; et al. Understanding the role of carbamate reactivity in fatty acid amide hydrolase inhibition by QM/MM mechanistic modelling. Chem. Commun. (Cambridge, U.K.) 2011, 9, 2517−2519. (19) Lodola, A.; Capoferri, L.; Rivara, S.; Tarzia, G.; Piomelli, D.; Mulholland, A.; Mor, M. Quantum mechanics/molecular mechanics modeling of fatty acid amide hydrolase reactivation distinguishes substrate from irreversible covalent inhibitors. J. Med. Chem. 2013, 56, 2500−2512. (20) Mileni, M.; Kamtekar, S.; Wood, D. C.; Benson, T. E.; Cravatt, B. F.; Stevens, R. C. Crystal structures of fatty acid amide hydrolase bound to the carbamate inhibitor URB597: discovery of a deacylating water molecule and insight into enzyme inactivation. J. Mol. Biol. 2010, 400, 743−754.

that the electrostatic term alone lacks the discriminating capability in this particular case because proper description of the FAAH inhibitory potency requires accounting for dispersive forces. Therefore, the combination of electrostatic and dispersion components seems to be the optimal measure of relative stabilities in biomolecular complexes. Related MP2 studies on the binding affinity of a hydrophobic p53−MDM2 oncoprotein complex79 or carbonic anhydrase inhibitors80 indicate the importance of dispersion forces. In the case of FAAH inhibition studied herein, several empirical scoring functions do not correlate with their experimental affinities. Overall, throughout the inhibitor design process, empirically derived results should be treated with caution as a negative outcome does not necessarily indicate a weak inhibitor. In terms of the highest percentage of correct results within all of the analyzed approaches, EMP2 is proven to be the most reliable but costly tool for FAAH inhibitor description. Because the computational cost associated with EMP2 calculations might be prohibitive, whenever large libraries of inhibitor-like compounds are concerned, the use of ab initio-derived atomic multipoles and dispersion atom−atom potential functions seems to be a reasonable alternative worthy of further testing.

■

ASSOCIATED CONTENT

* Supporting Information S

Structural features of the receptor−ligand models, performance of basis sets and approximate dispersion correction, CPU time requirements, and inhibitor binding characteristics. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was supported by the Wrocław Research Centre EIT+ under the project BIOMED “Biotechnologies and advanced medical technologies” (POIG 01.01.02-02-003/08) financed by the European Regional Development Fund Operational Programme Innovative Economy 1.1.2. We also thank the Wrocław University of Technology for support. The project was financed in part by a statutory activity subsidy from the Polish Ministry of Science and Higher Education for the faculty of Chemistry at the Wrocław University of Technology. Calculations were performed at the Wrocław Supercomputer and Networking Center (WCSS), the Poznań Supercomputer and Networking Center (PCSS), and the Interdisciplinary Center for Modeling (ICM) in Warsaw.

■

REFERENCES

(1) Lodola, A.; Mor, M.; Sirirak, J.; Mulholland, A. J. Insights into the mechanism and inhibition of fatty acid amide hydrolase from quantum mechanics/molecular mechanics (QM/MM) modelling. Biochem. Soc. Trans. 2009, 37, 363−367. (2) Hillard, C.; Weinlander, K.; Stuhr, K. Contributions of endocannabinoid signaling to psychiatric disorders in humans: gene and biochemical evidence. Neuroscience 2012, 204, 207−229. (3) Kerr, D.; Burke, N.; Ford, G.; Connor, T.; Harhen, B.; Egan, L.; Finn, D.; Roche, M. Pharmacological inhibition of endocannabinoid 14734

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

interaction to the prediction of binding affinity. J. Mol. Model. 2007, 13, 677−683. (43) Szarek, P.; Dyguda-Kazimierowicz, E.; Tachibana, A.; Sokalski, W. The physical nature of intermolecular interactions within cAMPdependent protein kinase active site: differential transition state stabilization in phosphoryl transfer reaction. J. Phys. Chem. B 2008, 112, 11819−11826. (44) Szefczyk, B.; Mulholland, A.; Ranaghan, K.; Sokalski, W. Differential transition-state stabilization in enzyme catalysis: quantum chemical analysis of interactions in the chorismate mutase reaction and prediction of the optimal catalytic field. J. Am. Chem. Soc. 2004, 126, 16148−16159. (45) Langner, K.; Kędzierski, P.; Sokalski, W.; Leszczyński, J. Physical nature of ethidium and proflavine interactions with nucleic acid bases in the interaction plane. J. Phys. Chem. B 2006, 110, 9720−9727. (46) Langner, K.; Janowski, T.; Góra, R.; Dziekoński, P.; Sokalski, W.; Pulay, P. The ethidium-UA/AU intercalation site: effect of model fragmentation and backbone charge state. J. Chem. Theory Comput. 2011, 7, 2600−2609. (47) Langner, K.; Beker, W.; Sokalski, W. Robust predictive power of the electrostatic term at shortened intermolecular distances. J. Phys. Chem. Lett. 2012, 3, 2785−2789. (48) Beker, W.; Langner, K.; Dyguda-Kazimierowicz, E.; Feliks, M.; Sokalski, W. Low cost prediction of relative stabilities of hydrogen bonded complexes from atomic multipole moments for overly short intermolecular distances. J. Comput. Chem. 2013, 34, 1797−1789. (49) Bracey, M. H.; Hanson, M. A.; Masuda, K. R.; Stevens, R. C.; Cravatt, B. F. Structural adaptations in a membrane enzyme that terminates endocannabinoid signaling. Science 2002, 298, 1793−1796. (50) Wu, Q.; Yang, W. T. Empirical correction to density functional theory for van der Waals interactions. J. Chem. Phys. 2002, 116, 515− 524. (51) Grimme, S. Accurate description of van der Waals complexes by density functional theory including empirical corrections. J. Comput. Chem. 2004, 25, 1463−1473. (52) Grimme, S. Semiempirical GGA-type density functional constructed with a long-range dispersion correction. J. Comput. Chem. 2006, 27, 1787−1799. (53) Korth, M. Error estimates for (semi-)empirical dispersion terms and large biomacromolecules. Org. Biomol. Chem. 2013, 11, 6515− 6519. (54) Pernal, K.; Podeszwa, R.; Patkowski, K.; Szalewicz, K. Dispersionless density functional theory. Phys. Rev. Lett. 2009, 103, 263201−263204. (55) Hoja, J.; Sax, A.; Szalewicz, K. Is electrostatics sufficient to describe hydrogen-bonding interactions? Chem.Eur. J. 2014, 20, 2292−2300. (56) McCammon, A. Theory of biomolecular recognition. Curr. Opin. Struct. Biol. 1998, 8, 245−249. (57) Gohlke, H.; Klebe, G. Approaches to the description and prediction of the binding affinity of small-molecule ligands to macromolecular receptors. Angew. Chem., Int. Ed. 2002, 41, 2644− 2676. (58) Gilson, M. K.; Zhou, H.-X. Calculation of protein-ligand binding affinities. Annu. Rev. Biophys. Biomol. Struct. 2007, 36, 21−42. (59) Olsson, T.; Williams, M.; Pitt, W.; Ladbury, J. The thermodynamics of protein-ligand interaction and solvation: insights for ligand design. J. Mol. Biol. 2008, 384, 1002−1017. (60) de M. Seabra, G.; Walker, R.; Roitberg, A. Are current semiempirical methods better than force fields? A study from the thermodynamics perspective. J. Phys. Chem. A 2009, 113, 11938− 11948. (61) Rao, L.; Zhang, I.; Guo, W.; Feng, L.; Meggers, E.; Xu, X. Nonfitting protein−ligand interaction scoring function based on firstprinciples theoretical chemistry methods: development and application on kinase inhibitors. J. Comput. Chem. 2013, 34, 1636−1646. (62) Kellogg, G.; Fornabaio, M.; Spyrakis, F.; Lodola, A.; Cozzini, P.; Mozzarelli, A.; Abraham, D. Getting it right. Modeling of pH, solvent

(21) Tarzia, G.; Duranti, A.; Gatti, G.; Piersanti, G.; Tontini, A.; Rivara, S.; Lodola, A.; Plazzi, P.; Mor, M.; Kathuria, S.; et al. Synthesis and structure-activity relationships of FAAH inhibitors: cyclohexylcarbamic acid biphenyl esters with chemical modulation at the proximal phenyl ring. ChemMedChem 2006, 1, 130−139. (22) Lodola, A.; Rivara, S.; Mor, M. Application of computational methods to the design of fatty acid amide hydrolase (FAAH) inhibitors based on a carbamic template structure. Adv. Protein Chem. Struct. Biol. 2011, 85, 1−26. (23) Sotriffer, C.; Stahl, M.; Klebe, G. Handbook of Chemoinformatics; Wiley-VC: Weinheim, Germany, 2003; pp 1732−1766. (24) Wolohan, P.; Reichert, D. CoMFA and docking study of novel estrogen receptor subtype selective ligands. J. Comput.-Aided Mol. Des. 2003, 17, 313−328. (25) Clark, R.; Norinder, U. Two personal perspectives on a key issue in contemporary 3D QSAR. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 108−113. (26) Stone, A. J. The Theory of Intermolecular Forces; Oxford University Press: Oxford, U.K., 1997. (27) Jeziorski, B.; Moszyński, R.; Szalewicz, K. Perturbation theory approach to intermolecular potential energy surfaces of van der Waals complexes. Chem. Rev. 1994, 94, 1887−1930. (28) Chalasinski, G.; Szcześniak, M. M. State of the art and challenges of the ab initio theory of intermolecular interactions. Chem. Rev. 2000, 100, 4227−4252. (29) Moszyński, R. Molecular Materials with Specific Interactions; Springer: Dordrecht, The Netherlands, 2007; pp 1−152. (30) Ř ezác,̌ J.; Riley, K.; Hobza, P. S66: A well-balanced database of benchmark interaction energies relevant to biomolecular structures. J. Chem. Theory Comput. 2011, 7, 2427−2438. (31) Frisch, M.; Trucks, G.; Schlegel, H.; Scuseria, G.; Robb, M.; Cheeseman, J.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G.; et al. Gaussian 09, revision D.01; Gaussian Inc.: Wallingford, CT, 2004. (32) Sokalski, W.; Kędzierski, P.; Grembecka, J. Ab initio study of physical nature of interactions between enzyme active site fragments in vacuo. Phys. Chem. Chem. Phys. 2001, 3, 657−663. (33) Boys, S.; Bernardi, F. The calculation of small molecular interactions by the differences of separate total energies. Some procedures with reduced errors. Mol. Phys. 1970, 19, 553−566. (34) Hariharan, P. C.; Pople, J. A. The influence of polarization functions on molecular orbital hydrogenation energies. Theor. Chim. Acta 1973, 28, 213−222. (35) Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; DeFrees, D. J.; Pople, J. A.; Gordon, M. S. Self-consistent molecular orbital methods. XXIII. A polarization-type basis set for second-row elements. J. Chem. Phys. 1982, 77, 3654−3665. (36) Rassolov, V. A.; Pople, J. A.; Ratner, M. A.; Windus, T. L. 631G* basis set for atoms K through Zn. J. Chem. Phys. 1998, 109, 1223−1229. (37) Góra, R.; Sokalski, W.; Leszczyński, J.; Pett, V. The nature of interactions in the ionic crystal of 3-pentenenitrile, 2-nitro-5-oxo, ion(−1), sodium. J. Phys. Chem. B 2005, 109, 2027−2033. (38) Schmidt, M.; Baldridge, K.; Boatz, J.; Elbert, S.; Gordon, M.; Jensen, J.; Koseki, S.; Matsunaga, N.; Nguyen, K.; Su, S.; et al. General atomic and molecular electronic structure system. J. Comput. Chem. 1993, 14, 1347−1363. (39) Woon, D. E.; Dunning, T. H., Jr. Gaussian basis sets for use in correlated molecular calculations. III. The atoms aluminum through argon. J. Chem. Phys. 1993, 98, 1358−1371. (40) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. Electron affinities of the first-row atoms revisited. Systematic basis sets and wave functions. J. Chem. Phys. 1992, 96, 6796−6806. (41) Dyguda, E.; Grembecka, J.; Sokalski, W.; Leszczyński, J. Origins of the activity of PAL and LAP enzyme inhibitors: toward ab initio binding affinity prediction. J. Am. Chem. Soc. 2005, 127, 1658−1659. (42) Grzywa, R.; Dyguda-Kazimierowicz, E.; Sieńczyk, M.; Feliks, M.; Sokalski, W.; Oleksyszyn, J. The molecular basis of urokinase inhibition: from the nonempirical analysis of intermolecular 14735

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736

The Journal of Physical Chemistry B

Article

and “nearly” everything else in virtual screening of biological targets. J. Mol. Graphics Modell. 2004, 22, 479−486. (63) Klebe, G. Virtual ligand screening: strategies, perspectives and limitations. Drug Discovery Today 2006, 11, 580−594. (64) Leach, A. R.; Shoichet, B. K.; Peishoff, C. E. Prediction of protein-ligand interactions. Docking and scoring: successes and gaps. J. Med. Chem. 2006, 49, 5851−5855. (65) Plewczynski, D.; Lazniewski, M.; Augustyniak, R.; Ginalski, K. Can we trust docking results? Evaluation of seven commonly used programs on PDBbind database. J. Comput. Chem. 2011, 32, 742−755. (66) Discovery Studio Modeling Environment, release 3.5; Accelrys Software Inc.: San Diego, 2012. (67) Tirado-Rives, J.; Jorgensen, W. Contribution of conformer focusing to the uncertainty in predicting free energies for proteinligand binding. J. Med. Chem. 2006, 49, 5880−5884. (68) Krammer, A.; Kirchhoff, P.; Jiang, X.; Venkatachalam, C.; Waldman, M. LigScore: a novel scoring function for predicting binding affinities. J. Mol. Graphics Modell. 2005, 23, 395−407. (69) Gehlhaar, D.; Verkhivker, G.; Rejto, P.; Sherman, C.; Fogel, D.; Fogel, L.; Freer, S. Molecular recognition of the inhibitor AG-1343 by HIV-1 protease: conformationally flexible docking by evolutionary programming. Chem. Biol. 1995, 2, 317−324. (70) Gehlhaar, D.; Bouzida, D.; Rejto, P. Rational Drug Design: Novel Methodology and Practical Applications; American Chemical Society: Washington, DC, 1999; Vol. 719, pp 292−311. (71) Jain, A. Scoring noncovalent protein-ligand interactions: a continuous differentiable function tuned to compute binding affinities. J. Comput.-Aided Mol. Des. 1996, 10, 427−440. (72) Muegge, I.; Martin, Y. A general and fast scoring function for protein-ligand interactions: a simplified potential approach. J. Med. Chem. 1999, 42, 791−804. (73) Muegge, I. PMF scoring revisited. J. Med. Chem. 2006, 49, 5895−5902. (74) Böhm, H. The development of a simple empirical scoring function to estimate the binding constant for a protein-ligand complex of known three-dimensional structure. J. Comput.-Aided Mol. Des. 1994, 8, 243−256. (75) Böhm, H. Prediction of binding constants of protein ligands: a fast method for the prioritization of hits obtained from the de novo design or 3D database search programs. J. Comput.-Aided Mol. Des. 1998, 12, 309−323. (76) Jorgensen, W. L. Efficient drug lead discovery and optimization. Acc. Chem. Res. 2009, 42, 724−733. (77) Xu, M.; Lill, M. A. Induced fit docking, and the use of QM/MM methods in docking. Drug Discovery Today: Technol. 2013, 10, e411− e418. (78) Langner, K.; Sokalski, W.; Leszczynski, J. J. Chem. Phys. 2007, 127, 111102−111104. (79) Ding, Y.; Mei, Y.; Zhang, J. Z. H. Quantum mechanical studies of residue-specific hydrophobic interactions in p53-MDM2 binding. J. Phys. Chem. B 2008, 112, 11396−11401. (80) Riley, K. E.; Cui, G.; Merz, K. M., Jr. An ab initio investigation of the interactions involving the aromatic group of the set of fluorinated N-(4-sulfamylbenzoyl)benzylamine inhibitors and human carbonic anhydrase II. J. Phys. Chem. B 2007, 111, 5700−5707.

14736

dx.doi.org/10.1021/jp5059287 | J. Phys. Chem. B 2014, 118, 14727−14736