A Fast Ab Initio Predictor Tool for Covalent ... - ACS Publications

Jun 27, 2019 - ABSTRACT: Thanks to their unique mode of action, covalent drugs represent an ... We tested this approach on acrylamides, one of the mos...
0 downloads 0 Views 1MB Size
Download PDF
Article Cite This: J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

pubs.acs.org/jcim

A Fast Ab Initio Predictor Tool for Covalent Reactivity Estimation of Acrylamides Ferruccio Palazzesi,* Marc A. Grundl, Alexander Pautsch, Alexander Weber, and Christofer S. Tautermann Medicinal Chemistry, Boehringer Ingelheim Pharma GmbH & Co. KG, Birkendorfer Strasse 65, 88397 Biberach an der Riss, Germany

Downloaded via GUILFORD COLG on July 17, 2019 at 13:57:17 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Thanks to their unique mode of action, covalent drugs represent an exceptional opportunity for drug design. After binding to a biologically relevant target system, covalent compounds form a reversible or irreversible covalent bond with a nucleophilic amino acid. Due to the inherently large binding energy of a covalent bond, covalent binders exhibit higher potencies and thus allow potentially lower drug dosages. However, a proper balancing of compound reactivity is key for the design of covalent binders, to achieve high levels of target inhibition while minimizing promiscuous covalent binding to nontarget proteins. In this work, we demonstrated the possibility to apply the electrophilicity index concept to estimate covalent compound reactivity. We tested this approach on acrylamides, one of the most prominent classes of covalent warheads. Our study clearly demonstrated that, for compounds with molecular weight (MW) below 250 Da, the electrophilicity index can be directly used to estimate compound reactivity. On the other hand, for leadlike molecules (MW > 250 Da) we implemented a new truncation algorithm that has to be applied before reactivity calculations. This algorithm can ensure the localization of HOMO/LUMO orbitals on the compound warhead and thus a correct estimation of its reactivity. Our results also indicate that caution should be used when employing the electrophilicity index to estimate the reactivity of nonterminal acrylamides. The nonparametric nature of this method and its reasonable computational cost make it a suitable tool to support covalent drug design.

■

INTRODUCTION In the past decade, the pharmaceutical industry has witnessed the successful resurgence of covalent drugs.1−5 A better understanding of their possible applicability in combination with their unique characteristics, such as higher potency, lower dosage, and longer duration of action, puts covalent drugs under a new light.1,4,5 For these reasons, the US Food and Drug Administration (FDA) already approved more than 50 covalent drugs for different therapeutic indications.1,6−8 Notably, three of them were present in the top-selling drugs list of the United States in 2009.1 As illustrated in Scheme 1, the standard mechanism of covalent inhibitor binding is a two-step process. While the first

step is influenced by noncovalent interactions between the compound and the target site, in the second step there is an interplay between proper compound positioning and the intrinsic reactivity of the compound reactive group, generally called the “warhead”, which finally leads to the formation of a covalent bond. One of the key aspects of the optimization process of covalent drugs is thus the modulation of reactivity of this warhead. This electrophilic group can react with specific nucleophilic residues, such as cysteine, lysine, tyrosine, or serine, and covalently bind the ligand to the target protein.6,9,10 Increasing the reactivity of the warhead results in faster formation of compound protein adducts, leading to improved potency. At the same time, selectivity and toxicity issues can occur due to nontarget specific, promiscuous binding of highly reactive compounds. Moreover, one must also consider an increased compound depletion by glutathione (GSH).1,5 Identifying and tuning the reactivity of suitable warheads is therefore fundamental for successfully optimizing covalent inhibitors.4 For this reason, researchers developed experimental assays to evaluate and classify the reactivity of covalent warheads.

Scheme 1. Standard, Two-Step, Covalent Mechanisma

a

E = enzyme, I = covalent inhibitor, E.I = noncovalent complex, E−I = covalent complex. While Ki depends only on noncovalent interactions, Kinact is also influenced by intrinsic compound reactivity. © XXXX American Chemical Society

Received: April 17, 2019 Published: June 27, 2019 A

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

suitable tool to calculate the covalent compound reactivity of huge virtual libraries in a time-efficient manner. Here, we decided to apply this approach to acrylamides, a prominent warhead class that is already well-established in the market.1,37 Moreover, this warhead type can span a wide range of reactivity,26 making it a perfect model for our study. Due to the fact that the electrophilicity index is based on highest occupied (HOMO) and lowest unoccupied (LUMO) molecular orbitals energies to estimate compound reactivity, these frontier orbitals must be localized on the warhead. This means that some precautions have to be taken to ensure a correct localization of these orbitals. While for fragment-size compounds (MW < 250 Da38) no modifications are needed, we solved this issue for larger compounds through a truncation algorithm. One can argue that compound reactivity can also be influenced by the presence of specific interactions and/or steric clashes with the target system, such as the GSH peptide or the protein binding site. Here, in order to keep our approach more general, we calculate compound reactivity only on the basis of its electrophilicity, without including such target-specific interaction effects.

Adduct formation reaction with GSH peptide is the standard way to evaluate the reactivity of compounds toward cysteines.11−15 In this assay, an excess of GSH is dissolved in aqueous phosphate buffer solution at physiological pH. After the addition of an electrophilic compound, the kinetics of adduct formation can be measured to assess compound reactivity. Recently, a similar probe (N-α-acetyl-L-lysine, hereafter LYS) has been reported to assess compound reactivity toward lysine.16,17 While such reactivity assays are hence useful, the data generated in all these measurements highly depend on several factors, such as experimental conditions, compound degradation in buffer, and reversibility of covalent bond formation. This makes a direct comparison between compounds, especially when tested in different laboratories, extremely difficult, as also recently highlighted by Baillie et al.5 Moreover, one has to consider that the timeand capacity-consuming preparation of compounds is always required before reactivity can be assessed experimentally. In this context, a more reliable and general method to predict covalent warhead reactivity would be beneficial, supporting the design of new compounds. In order to overcome all the aforementioned experimental issues, we decided to develop a computational approach to predict compound reactivity. Several different computational methods able to predict the reactivity of electrophilic compounds can be found in the literature. A standard tool nowadays used in computer-aided drug design is the quantitative structure−activity relationship (QSAR) approach. One of the simplest QSAR models that can be used to estimate experimental GSH data is based on Hammett’s equation.18 However, this model only allows studying the effects of aromatic ring substitutions.19 Other more advanced QSAR approaches are present in the literature; unfortunately, very few of them take into consideration covalent reactivity and related toxicity.20−23 Moreover, the scarce availability of experimental GSH half-life (t1/2) data and their aforementioned experimental variation and uncertainty5,15 does not allow building QSAR models with broad applicability. Recently, Hughes et al. trained a neural network with more than 1000 different molecules known to be reactive toward GSH.24 However, this model did not incorporate any quantitative experimental GSH data, limiting model performances to predict compound reactivity. Accurate estimation of compound reactivity with respect to model cysteines can be obtained from transition-state calculations and/or adduct formation energy estimations by quantum chemistry.14,15,25−28 However, all these methodologies require a priori knowledge of reaction pathways or the application of algorithms to locate transition states. Moreover, their implementation in a highthroughput fashion is not straightforward, mainly due to their computational complexity. To overcome all these limitations, we decided to investigate the application of the “quantum chemical electrophilicity index” (ω) developed by Parr et al.29 for covalent warhead reactivity. This ab initio approach, or its atom-localized version,30,31 has successfully been applied to calculate compound electrophilicity and toxicity, alone or in combination with other descriptors.20,22,24,32−36 The electrophilicity index can be easily calculated by density functional theory (DFT) calculations. Moreover, this approach is very general and thus can be applied to different types of compound classes, in contrast to the decision-tree-based approach recently suggested by Lonsdale et al.26 All these properties make it a

■

MATERIALS AND METHODS Quantum Mechanical (QM) Calculations. The electrophilicity index approach is based on the electronic energies of LUMO (ELUMO) and HOMO (EHOMO) molecular orbitals:29,39 ω=

μ2 2η

(1)

where μ≈

E HOMO + E LUMO 2

η ≈ E LUMO − E HOMO

(electronic chemical potential)

(chemical hardness)

For a thorough explanation of this approach, please refer to refs 40 and 41. Starting from the 3D compound configuration generated with CORINA,42−44 ω is calculated using eq 1) after geometry optimization at the B3LYP/6-31+G(d,p) level of theory. To take into account solvent effects, we considered implicit water using CPCM.45,46 For conformational search analysis, we used OMEGA conformational search47,48 with an energy window of 10 kcal/mol and a maximum of 200 conformations per compound. We performed a QM geometry optimization for all the compound conformations and then used the lowest QM energy conformation to calculate ω. All the DFT calculations have been performed with Gaussian 09.49 To assess the correlation of the calculated reactivity with respect to experimental GSH data, we calculated the squared correlation coefficient (R2) and Spearman coefficient (S). While the former one measures the quality of a linear correlation between data, the latter one measures the rank correlation. From both correlation analyses we excluded compounds for which experimental GSH data are reported with operator values (> or < ). Experimental GSH and LYS adduct formation data were taken from refs 15, 16, 19, 26, and 50. Compound structures and experimental data are reported in Table S4 of Supporting Information (SI). B

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling Truncation Algorithm. In order to correctly localize HOMO/LUMO orbitals on the compound warhead, and thus to estimate covalent compound reactivity, here we introduced a new truncation procedure. This algorithm starts by including all atoms up to three bonds away from the nitrogen and the βcarbon atom of the acrylamide moiety. In the case in which one of these atoms is part of a functional group, we include the whole functional group. For the full list of considered functional groups, see Table S1 (SI). We must point out that ring systems are also considered as functional groups. In this way, if one atom of the three-bond shell belongs to a ring (aromatic or not), we include the whole ring system to the final truncated compound. In the case in which ring substituents are present, we also include them according to the aforementioned functional groups and rings dictionary. A schematic flowchart of this algorithm can be found in Scheme 2. This truncation algorithm has been implemented using RDKit functionalities.51

sources cannot be pooled together due to the high variability of the GSH adduct formation assay, as mentioned before. Fragmentlike Compounds. In Figure 1a we report the electrophilicity index calculated for compounds 1−14, 29, and 30 (see Table S4, SI).26 This ab initio quantity correlates well with GSH experimental data: R2 of 0.78 and Spearman correlation of 0.96. A similar correlation has been recently observed for the same compound set by Lonsdale et al.26 using the DFT activation energy barrier and by Ward et al.50 with the Hammett approach. In Table 1 we summarize different methods’ performances. Dahal et al. established an assay to estimate compound reactivity to lysine using N-α-acetyl-L-lysine as model peptide.16,17 In Figure 1b we report the electrophilicity index calculated for compounds 109−143, excluding compound 127 because it does not contain an acrylamide warhead. This result demonstrates that ω can also be used to correctly estimate acrylamide reactivity with respect to lysine (R2 = 0.61 and S = 0.74). We must point out that outliers are present in both analyses of Figure 1. In particular, the presence of substitutions on the reactive double bond of acrylamide (see squares in Figure 1) leads to a deviation from linear correlation. For these cases and in accordance with literature data,15,16,52 experimental GSH reactivity is also clearly influenced by steric hindrance, which is not considered in the electrophilicity index approach.53 We obtained a similar trend using ω also with other series of acrylamide fragments for which experimental GSH data are available (see Table S4, SI). Flanagan et al.15 calculated the activation free energies to estimate covalent reactivity, obtaining a R2 with respect to the experimental data of 0.915. As reported in Table 1 and Figure S1 (SI), our calculations show a R2 of 0.68, with a Spearman correlation of 0.91. Cee et al.19 studied the ability of Hammett parametrization to reproduce ortho, meta, and para substitutions effects on GSH reactivity. The reported R2 correlation with respect to the experimental data is equal to 0.80 (ortho), 0.83 (meta), and 0.98 (para). The authors of this paper compared these results with correlation obtained with activation barrier energies calculated at the DFT level. The reported R2 values for activation energy and activation free energy are 0.83 and 0.75, respectively. The electrophilicity index gave us a R2 coefficient of 0.64 and a Spearman coefficient of 0.88 considering all the ortho, meta, and para substituents together (see Table 1). As can be seen from Figure S2 (SI), compounds possessing a nitro group as a ring substituent (38−40) lead to high calculated reactivity values. Upon removal of these values, correlation with experimental GSH data further improved (R2 = 0.79, S = 0.85). All these results suggest that the electrophilicity index can achieve a good accuracy with respect to the Hammett method and a lower, but still remarkable, agreement with transition state (TS) energies. Furthermore, it is more general than the Hammett model, and a computationally less time demanding method compared to activation barrier energies. The calculation of the TS and the required vibrational frequencies can be in fact 10 times more demanding than doing the energy minimization run needed for the electrophilicity index calculation.54 To further check the robustness of our results, we calculated ω for compounds 1−14, 29, and 30 using different DFT functionals and basis sets. The agreement with experimental data did not improve with respect to that obtained with

Scheme 2. Truncation Algorithm Developed To Study Covalent Compound Reactivity and the Effects of This on Compound 19 of Ref 26a

a

For the full list of functional groups and ring systems considered, see Table S1 (SI).

■

RESULTS AND DISCUSSION As introduced before, in this work we make use of ω to calculate covalent compound reactivity. This approach is based on frontier orbitals, HOMO and LUMO, and thus, these must be located on the covalent warhead in order to correctly estimate compound reactivity. Therefore, after collecting compounds and associated experimental GSH data from the literature (see Materials and Methods), we decided to split them into “fragmentlike” and “leadlike” molecules according to their size. All compounds with MW below 250 Da38 have been assigned to the first group and all the others to the second one. We started our investigation by calculating ω for compounds belonging to the first family. We would like to highlight the fact that compounds with experimental data coming from different C

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

Figure 1. Correlation between calculated reactivity (ω) and experimental data. Correlation with respect to (a) GSH data and (b) LYS data. Circles indicate unsubstituted acrylamide compounds, while squares highlight substituted ones. The color code is based on experimental data: green, high reactivity (t1/2 < 1 h); yellow, moderate reactivity (1 h < t1/2 < 10 h); orange, low reactivity (10 h < t1/2 < 100 h); and red, not reactive at all. Experimental GSH/LYS adduct formation t1/2 axis is in log units.

Table 1. Reactivity Trends As Described by the Calculated Electrophilicity Index and Other Methods this work data set

compounds

Lonsdale et al.26

1−10, 29, 30

Ward et al.50

1−10, 29, 30

Hammett R2 S R2 S

Flanagan et al.15

88−93, 96, 102−108

Cee et al.19

31−64c

R2 S R2

S

Dahal et al.16

109−143

ΔGa barrier

a

ΔEa barrier

a

0.75 0.86 0.92 0.91 0.94 0.92

0.80 0.83 0.98 0.89 0.91 0.98

(m) (i) (m) (i)

results 0.78 0.96 0.78

data Figure 1a Figure 1a

0.96 0.86 (0.915)b 0.96 0.75

(o) (m) (p) (o) (m) (p)

0.83

0.62 (0.68)b 0.89 (0.91)b 0.64 (0.79)d

Figure S1

Figure S2 0.87

R2 S

0.91

0.88 (0.85)d

0.61 0.74

Figure 1b

ΔGa and ΔEa are the activation free energy and activation energy, respectively. bNote that R2 and S differ from the aforementioned ones (in parentheses) because values for the ΔGa barrier of compounds 97 and 98 are not reported in the original paper. Thus we compared R2 and S using compounds 88−93, 96, and 102−108. cCompound 31 is not considered in the original paper because it does not possess any aromatic ring substitutions. dCorrelations obtained after removing compounds containing nitro groups as ring substituents (see the text). a

B3LYP/6-31+G(d,p) (see Table S2, SI). To check whether the electrophilicity estimation is dependent on the starting compound conformations, we performed a thorough conformational search using compounds 1−14, 29, and 30 (see Materials and Methods for further details). Calculating ω from the lowest QM energy conformation, we obtained a similar agreement with experimental GSH data compared to that where only a single conformation has been used (see Figure S3 and Table S3, SI). Leadlike Compounds. We then decided to extend the electrophilicity index calculations also to larger compounds (MW > 250 Da). Different from fragmentlike compounds, multiple effects can have an influence on the final warhead reactivity in these molecules. For example, the presence of several reactive groups in the compound, far or close to the warhead, and/or the combination of electron-withdrawing and -donating substituents can affect compound reactivity. Before calculating ω for such kind of compounds we decided to analyze the localization of frontier orbitals. As can be seen for compound 19 in Figure 2a, the LUMO orbital is localized on

the aromatic part of the molecule, far from the acrylamide moiety [see Figure S4 (SI) for the HOMO orbital]. Lonsdale et al. observed a similar delocalization when trying to use the LUMO energy to estimate covalent compound reactivity.26 As seen in eq 1, the electrophilicity index theory is based on the HOMO and LUMO energies, and thus, the localization of these orbitals far apart from the covalent warhead will lead to wrong reactivity estimations, resulting in poor correlation with experimental GSH data. This is what we observed in Figure 2c, in which we report ω calculated for leadlike compounds 15− 28 (R2 = 0.30, Spearman correlation = 0.48). Here, we did not consider compound 20 because the basis set used in this study cannot be used for iodine, while compound 25 has been removed because it contains two different warheads (acrylamide and vinylsulfonamide). However, it is reasonable to suppose that not all the atoms in the molecule will influence compound warhead reactivity.26,55 This assumption can be easily verified by comparing available experimental data with compound structure. As can be seen in Figure S5 (SI), all these compounds (65−68, 73, D

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

Figure 2. The upper panels show the graphical representation of the LUMO orbital for compound 19. (a) LUMO localization considering the whole molecule. (b) LUMO position after molecule truncation. The bottom panels show the electrophilicity index calculated for leadlike compounds 15−28 plotted against experimental GSH data. Results considering (c) the whole molecule and (d) after applying truncation rules described in Scheme 2. Circles indicate unsubstituted acrylamide compounds, while squares highlight substituted ones. The color code is based on experimental data: green, high reactivity (t1/2 < 1 h); yellow, moderate reactivity (1 h < t1/2 < 10 h); orange, low reactivity (10 h < t1/2 < 100 h); and red, not reactive at all. Experimental GSH adduct formation t1/2 axis is in log units.

■

CONCLUSION In this work, we extend the application of the electrophilicity index idea to estimate the absolute covalent warhead reactivity. We apply this to acrylamides, a prominent warhead class already present on the market with covalent drugs such as afatinib, ibrutinib, osimertinib, and others. The results reported here demonstrated good correlations between the available experimental GSH reactivity data and the calculated electrophilicity index. A similar outcome is also observed with respect to experimental reactivity data using lysine as the target amino acid. While for compounds with MW below 250 Da the application of such a computational approach is rather straightforward, caution must be taken for leadlike molecules. A localization of HOMO/LUMO orbitals far from the compound warhead can in fact preclude the correct estimation of compound reactivity. For this reason, we introduced here a new algorithm to truncate molecules that allows one to correctly localize these orbitals on the warhead moiety. After the application of this algorithm, the electrophilicity index correctly estimated covalent warhead reactivity also for such bigger compounds. We must also admit that the electrophilicity index approach could not always correctly predict the reactivity of nonterminal acrylamide compounds. For this class of compounds, steric clashes and/or specific interactions with the target system could indeed strongly influence covalent warhead reactivity, and this effect could not be accurately predicted by our approach. The possibility to obtain information on compound reactivity in such a fast manner will open new opportunities for covalent drug design, in order to obtain safer and more potent covalent compounds.

75, 80) possess a methoxyaniline-substituted phenyl ring attached to the acrylamide nitrogen atom and share similar experimental GSH data. For this reason, we decided to truncate compounds according to some simple chemical rules explained in Materials and Methods and summarized in Scheme 2. After the application of this truncation algorithm, we analyzed again the LUMO and HOMO localization on the compounds. As can be seen from Figures 2b and S4 (SI), now the frontier orbitals are occupying the correct position to calculate compound reactivity. We thus calculated ω of leadlike molecules 15−28,26 after applying this truncation algorithm. Once again we excluded compound 25 for the presence of two different covalent warheads (acrylamide and vinylsulfonamide). As can be seen from Figure 2d, the electrophilicty index calculated now correlates well with experimental GSH data (R2 = 0.91, S = 0.84). In this case, ω also correctly predicted the reactivity of a nonterminal acrylamide such as compound 24. In Figure S6 (SI), we report another set of molecules taken from ref 50 for which compound reactivity can be better estimated after truncation (R2 from 0.23 to 0.52 and S from 0.54 to 0.63). Some of these compounds possess a MW above 250 Da also after the application of the truncation algorithm (see for example 69 and 81−86). However, a repeated application of this procedure will give the same compound structure. We can conclude that the electrophilicity index in combination with the truncation algorithm can be considered an accurate and efficient method to estimate the absolute acrylamide warhead reactivity of leadlike molecules. E

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

■

(9) Shannon, D. A.; Weerapana, E. Covalent Protein Modification: the Current Landscape of Residue-Specific Electrophiles. Curr. Opin. Chem. Biol. 2015, 24, 18−26. (10) Mukherjee, H.; Grimster, N. P. Beyond Cysteine: Recent Developments in the Area of Targeted Covalent Inhibition. Curr. Opin. Chem. Biol. 2018, 44, 30−38. (11) Bohme, A.; Thaens, D.; Paschke, A.; Schuurmann, G. Kinetic Glutathione Chemoassay to Quantify Thiol Reactivity of Organic Electrophiles–Application to Alpha, Beta-Unsaturated Ketones, Acrylates, and Propiolates. Chem. Res. Toxicol. 2009, 22, 742−750. (12) Schwöbel, J. A.; Koleva, Y. K.; Enoch, S. J.; Bajot, F.; Hewitt, M.; Madden, J. C.; Roberts, D. W.; Schultz, T. W.; Cronin, M. T. Measurement and Estimation of Electrophilic Reactivity for Predictive Toxicology. Chem. Rev. 2011, 111, 2562−2596. (13) MacFaul, P. A.; Morley, A. D.; Crawford, J. J. A Simple in Vitro Assay for Assessing the Reactivity of Nitrile Containing Compounds. Bioorg. Med. Chem. Lett. 2009, 19, 1136−1138. (14) Berteotti, A.; Vacondio, F.; Lodola, A.; Bassi, M.; Silva, C.; Mor, M.; Cavalli, A. Predicting the Reactivity of Nitrile-Carrying Compounds with Cysteine: a Combined Computational and Experimental Study. ACS Med. Chem. Lett. 2014, 5, 501−505. (15) Flanagan, M. E.; Abramite, J. A.; Anderson, D. P.; Aulabaugh, A.; Dahal, U. P.; Gilbert, A. M.; Li, C.; Montgomery, J.; Oppenheimer, S. R.; Ryder, T.; et al. Chemical and Computational Methods for the Characterization of Covalent Reactive Groups for the Prospective Design of Irreversible Inhibitors. J. Med. Chem. 2014, 57, 10072−10079. (16) Dahal, U. P.; Gilbert, A. M.; Obach, R. S.; Flanagan, M. E.; Chen, J. M.; Garcia-Irizarry, C.; Starr, J. T.; Schuff, B.; Uccello, D. P.; Young, J. A. Intrinsic Reactivity Profile of Electrophilic Moieties to Guide Covalent Drug Design: N-α-acetyl-L-lysine as an Amine Nucleophile. MedChemComm 2016, 7, 864−872. (17) Mukherjee, H.; Debreczeni, J.; Breed, J.; Tentarelli, S.; Aquila, B.; Dowling, J.; Whitty, A.; Grimster, N. A Study of the Reactivity of S (VI)−F Containing Warheads with Nucleophilic Amino-Acid Side Chains Under Physiological Conditions. Org. Biomol. Chem. 2017, 15, 9685−9695. (18) Hansch, C.; Leo, A.; Taft, R. A Survey of Hammett Substituent Constants and Resonance and Field Parameters. Chem. Rev. 1991, 91, 165−195. (19) Cee, V. J.; Volak, L. P.; Chen, Y.; Bartberger, M. D.; Tegley, C.; Arvedson, T.; McCarter, J.; Tasker, A. S.; Fotsch, C. Systematic Study of the Glutathione (GSH) Reactivity of N-Arylacrylamides: 1. Effects of Aryl Substitution. J. Med. Chem. 2015, 58, 9171−9178. (20) Schwöbel, J. A.; Wondrousch, D.; Koleva, Y. K.; Madden, J. C.; Cronin, M. T.; Schüürmann, G. Prediction of Michael-Type Acceptor Reactivity Toward Glutathione. Chem. Res. Toxicol. 2010, 23, 1576− 1585. (21) Shokhen, M.; Traube, T.; Vijayakumar, S.; Hirsch, M.; Uritsky, N.; Albeck, A. Differentiating Serine and Cysteine Protease Mechanisms by new Covalent QSAR Descriptors. ChemBioChem 2011, 12, 1023−1026. (22) Smith, J. M.; Rowley, C. N. Automated Computational Screening of the Thiol Reactivity of Substituted Alkenes. J. Comput.Aided Mol. Des. 2015, 29, 725−735. (23) Chan, K.; Jensen, N.; O’Brien, P. J. Structure−Activity Relationships for Thiol Reactivity and Rat or Human Hepatocyte Toxicity Induced by Substituted p Benzoquinone Compounds. J. Appl. Toxicol. 2008, 28, 608−620. (24) Hughes, T. B.; Miller, G. P.; Swamidass, S. J. Site of Reactivity Models Predict Molecular Reactivity of Diverse Chemicals with Glutathione. Chem. Res. Toxicol. 2015, 28, 797−809. (25) Schwöbel, J.; Madden, J.; Cronin, M. Examination of Michael Addition Reactivity Towards Glutathione by Transition-State Calculations. SAR QSAR Environ. Res. 2010, 21, 693−710. (26) Lonsdale, R.; Burgess, J.; Colclough, N.; Davies, N. L.; Lenz, E. M.; Orton, A. L.; Ward, R. A. Expanding the Armory: Predicting and Tuning Covalent Warhead Reactivity. J. Chem. Inf. Model. 2017, 57, 3124−3137.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jcim.9b00316. Functional groups dictionary considered in the truncation algorithm, the electrophilicity index results in comparison with other in silico methods to predict covalent warhead reactivity, influence on the electrophilicity index of different basis set/functional combinations and conformational search analysis, HOMO orbital location before and after truncation, comparison between compound structures and experimental GSH data obtained from the literature, the electrophilicity index calculated for another set of leadlike compounds, and list of compound structures employed in this work and associated experimental GSH and LYS adduct formation data. (PDF)

■

Simplified molecular-input line-entry system (SMILES) notations for the whole molecules and fragments used (XLSX)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Ferruccio Palazzesi: 0000-0003-1534-4781 Christofer S. Tautermann: 0000-0002-6935-6940 Author Contributions

All authors have equally contributed to this work. All authors have given approval to the final version of the manuscript. Notes

The authors declare no competing financial interest.

■

ABBREVIATIONS GSH, glutathione; MW, molecular weight; LYS, N-α-acetyl-Llysine; DFT, density functional theory; QM, quantum mechanical; TS, transition state

■

REFERENCES

(1) Singh, J.; Petter, R. C.; Baillie, T. A.; Whitty, A. The Resurgence of Covalent Drugs. Nat. Rev. Drug Discovery 2011, 10, 307−317. (2) Johnson, D. S.; Weerapana, E.; Cravatt, B. F. Strategies for Discovering and Derisking Covalent, Irreversible Enzyme Inhibitors. Future Med. Chem. 2010, 2, 949−964. (3) Strelow, J. M. A Perspective on the Kinetics of Covalent and Irreversible Inhibition. SLAS Discov 2017, 22, 3−20. (4) De Cesco, S.; Kurian, J.; Dufresne, C.; Mittermaier, A. K.; Moitessier, N. Covalent Inhibitors Design and Discovery. Eur. J. Med. Chem. 2017, 138, 96−114. (5) Baillie, T. A. Targeted Covalent Inhibitors for Drug Design. Angew. Chem., Int. Ed. 2016, 55, 13408−13421. (6) Gehringer, M.; Laufer, S. A. Emerging and Re-Emerging Warheads for Targeted Covalent Inhibitors: Applications in Medicinal Chemistry and Chemical Biology. J. Med. Chem. 2019, 62, 5673. (7) Bauer, R. A. Covalent Inhibitors in Drug Discovery: From Accidental Discoveries to Avoided Liabilities and Designed Therapies. Drug Discovery Today 2015, 20, 1061−1073. (8) Vasudevan, A.; Argiriadi, M. A.; Baranczak, A.; Friedman, M. M.; Gavrilyuk, J.; Hobson, A. D.; Hulce, J. J.; Osman, S.; Wilson, N. S. Covalent Binders in Drug Discovery. Prog. Med. Chem. 2019, 58, 1− 62. F

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling (27) Mulliner, D.; Wondrousch, D.; Schüürmann, G. Predicting Michael-Acceptor Reactivity and Toxicity through Quantum Chemical Transition-State Calculations. Org. Biomol. Chem. 2011, 9, 8400− 8412. (28) Oballa, R. M.; Truchon, J.-F.; Bayly, C. I.; Chauret, N.; Day, S.; Crane, S.; Berthelette, C. A Generally Applicable Method for Assessing the Electrophilicity and Reactivity of Diverse NitrileContaining Compounds. Bioorg. Med. Chem. Lett. 2007, 17, 998− 1002. (29) Parr, R. G.; Szentpaly, L. v.; Liu, S. Electrophilicity Index. J. Am. Chem. Soc. 1999, 121, 1922−1924. (30) Parr, R. G.; Yang, W. Density Functional Approach to the Frontier-Electron Theory of Chemical Reactivity. J. Am. Chem. Soc. 1984, 106, 4049−4050. (31) Yang, W.; Mortier, W. J. The use of Global and Local Molecular Parameters for the Analysis of the Gas-Phase Basicity of Amines. J. Am. Chem. Soc. 1986, 108, 5708−5711. (32) Pereira, F.; Latino, D. A.; Aires-De-Sousa, J. Estimation of Mayr Electrophilicity with a Quantitative Structure−Property Relationship Approach using Empirical and DFT Descriptors. J. Org. Chem. 2011, 76, 9312−9319. (33) Kostal, J.; Voutchkova-Kostal, A. CADRE-SS, an in Silico Tool for Predicting Skin Sensitization Potential Based on Modeling of Molecular Interactions. Chem. Res. Toxicol. 2016, 29, 58−64. (34) Enoch, S. J.; Cronin, M. T. D.; Schultz, T. W.; Madden, J. C. Quantitative and Mechanistic Read Across for Predicting the Skin Sensitization Potential of Alkenes Acting via Michael Addition. Chem. Res. Toxicol. 2008, 21, 513−520. (35) Wondrousch, D.; Böhme, A.; Thaens, D.; Ost, N.; Schüürmann, G. Local Electrophilicity Predicts the Toxicity-Relevant Reactivity of Michael Acceptors. J. Phys. Chem. Lett. 2010, 1, 1605−1610. (36) Schwöbel, J. A.; Koleva, Y. K.; Enoch, S. J.; Bajot, F.; Hewitt, M.; Madden, J. C.; Roberts, D. W.; Schultz, T. W.; Cronin, M. T. Measurement and estimation of electrophilic reactivity for predictive toxicology. Chem. Rev. 2011, 111, 2562−2596. (37) Zhao, Z.; Bourne, P. E. Progress with Covalent Small-Molecule Kinase Inhibitors. Drug Discovery Today 2018, 23, 727−735. (38) Erlanson, D. A. Introduction to Fragment-Based Drug Discovery. In Fragment-based drug discovery and X-ray crystallography; Springer, 2011; pp 1−32. (39) Koopmans, T. Ü ber die Zuordnung von Wellenfunktionen und Eigenwerten zu den Einzelnen Elektronen eines Atoms. Physica 1934, 1, 104−113. (40) Chattaraj, P. K.; Giri, S.; Duley, S. Update 2 of: Electrophilicity Index. Chem. Rev. 2011, 111, PR43−PR75. (41) Domingo, L.; Ríos-Gutiérrez, M.; Pérez, P. Applications of the Conceptual Density Functional Theory Indices to Organic Chemistry Reactivity. Molecules 2016, 21, 748. (42) Sadowski, J.; Gasteiger, J.; Klebe, G. Comparison of Automatic Three-Dimensional Model Builders using 639 X-ray Structures. J. Chem. Inf. Model. 1994, 34, 1000−1008. (43) Schwab, C. H. Conformations and 3D Pharmacophore Searching. Drug. Discovery Today: Technol. 2010, 7, e245−e253. (44) 3D Structure Generator CORINA Classic; MN-AM: Nuremberg, Germany (www.mn-am.com). (45) Barone, V.; Cossi, M. Quantum Calculation of Molecular Energies and Energy Gradients in Solution by a Conductor Solvent Model. J. Phys. Chem. A 1998, 102, 1995−2001. (46) Cossi, M.; Rega, N.; Scalmani, G.; Barone, V. Energies, Structures, and Electronic Properties of Molecules in Solution with the C PCM Solvation Model. J. Comput. Chem. 2003, 24, 669−681. (47) Hawkins, P. C.; Skillman, A. G.; Warren, G. L.; Ellingson, B. A.; Stahl, M. T. Conformer Generation with OMEGA: Algorithm and Validation Using High Quality Structures from the Protein Databank and Cambridge Structural Database. J. Chem. Inf. Model. 2010, 50, 572−584. (48) Nm, S. F. OMEGA 3.1.0.3; Open Eye Scientific Software (http://www.eyesopen.com).

(49) Frisch, M.; Trucks, G.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G.; et al. Gaussian 09; Gaussian, Inc.: Wallingford, CT, 2009. (50) Ward, R. A.; Anderton, M. J.; Ashton, S.; Bethel, P. A.; Box, M.; Butterworth, S.; Colclough, N.; Chorley, C. G.; Chuaqui, C.; Cross, D. A.; et al. Structure-and Reactivity-Based Development of Covalent Inhibitors of the Activating and Gatekeeper Mutant Forms of the Epidermal Growth Factor Receptor (EGFR). J. Med. Chem. 2013, 56, 7025−7048. (51) Landrum, G. RDKit: Cheminformatics and Machine Learning Software; 2016 (http://www.rdkit.org). (52) Jackson, P. A.; Widen, J. C.; Harki, D. A.; Brummond, K. M. Covalent Modifiers: a Chemical Perspective on the Reactivity of α, βUnsaturated Carbonyls with Thiols via Hetero-Michael Addition Reactions. J. Med. Chem. 2017, 60, 839−885. (53) Allgäuer, D. S.; Jangra, H.; Asahara, H.; Li, Z.; Chen, Q.; Zipse, H.; Ofial, A. R.; Mayr, H. Quantification and Theoretical Analysis of the Electrophilicities of Michael Acceptors. J. Am. Chem. Soc. 2017, 139, 13318−13329. (54) Jensen, F. Introduction to Computational Chemistry; John Wiley & Sons, 1999. (55) Ebbrell, D. J.; Madden, J. C.; Cronin, M. T.; Schultz, T. W.; Enoch, S. J. Development of a Fragment-Based in Silico Profiler for Michael Addition Thiol Reactivity. Chem. Res. Toxicol. 2016, 29, 1073−1081.

G

DOI: 10.1021/acs.jcim.9b00316 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX