Computing pKa Values in Different Solvents by ... - ACS Publications

Article pubs.acs.org/JCTC

Computing pKa Values in Different Solvents by Electrostatic Transformation Emanuele Rossini,† Roland R. Netz,‡ and Ernst-Walter Knapp*,† †

Institute of Chemistry and Biochemistry, Freie Universität Berlin, Fabeckstrasse 36a, D-14195 Berlin, Germany Department of Physics, Freie Universität Berlin, Arnimallee 14, D-14195 Berlin, Germany

‡

S Supporting Information *

ABSTRACT: We introduce a method that requires only moderate computational effort to compute pKa values of small molecules in different solvents with an average accuracy of better than 0.7 pH units. With a known pKa value in one solvent, the electrostatic transform method computes the pKa value in any other solvent if the proton solvation energy is known in both considered solvents. To apply the electrostatic transform method to a molecule, the electrostatic solvation energies of the protonated and deprotonated molecular species are computed in the two considered solvents using a dielectric continuum to describe the solvent. This is demonstrated for 30 molecules belonging to 10 different molecular families by considering 77 measured pKa values in 4 different solvents: water, acetonitrile, dimethyl sulfoxide, and methanol. The electrostatic transform method can be applied to any other solvent if the proton solvation energy is known. It is exclusively based on physicochemical principles, not using any empirical fetch factors or explicit solvent molecules, to obtain agreement with measured pKa values and is therefore ready to be generalized to other solute molecules and solvents. From the computed pKa values, we obtained relative proton solvation energies, which agree very well with the proton solvation energies computed recently by ab initio methods, and used these energies in the present study.

1. INTRODUCTION The pH-dependent protonation state of titratable molecules in different solvents is given by the pKa value. It is a physicochemical key quantity to characterize and predict the behavior of such molecules.1−3 The pKa governs the solubility, volume of distribution and clearance of drug molecules,4−8 and the thermodynamics of metabolic reactions.9 Drugs usually need to penetrate cell membranes to become active, which they can do efficiently only in the uncharged state.7 The charge state of a drug depends on its pKa value and environment. Cell membranes are less polar than water. Hence, knowing the pKa value of a drug in nonpolar organic solvents can be very helpful. Having such information of a molecule at hand before it is considered seriously as a potential drug and synthesized can save a lot of tedious work of organic chemists and pharmacists. This knowledge is also useful for physicochemists when choosing the appropriate combination of a titratable molecule and a solvent to obtain a defined protonation state. Appropriate and inexpensive computational methods can help to design a desired pKa value for a titratable molecule belonging to a specific family. There are many methods for computing or predicting pKa values of small molecules solvated in water. The methods vary from mainly empirical over mixed to essentially ab initio approaches. The latter approaches usually use a combination of quantum chemistry (QC), molecular dynamics (MD), and © 2016 American Chemical Society

electrostatics as reviewed in refs 10−20. The empirical methods use machine learning tools and employ molecular descriptors, which can be based on molecular structure information alone21−23 or on a combination of structure and semiempirical QC24−30 or on other quantities such as group philicities.31 In the mixed approaches, an approximate pKa value is computed by a combination of QC and electrostatic solvation energies (ESEs), which enters a linear regression function, whose parameters are adjusted and may depend on molecular families.32−39 There are also promising approaches that essentially use classical MD simulations.40−47 There are different ab initio methods to compute the pKa value of a titratable group in water using an explicit representation of the solvent. Some describe the whole molecular system quantum mechanically48−50 or use for instance the Car−Parrinello approach.51 Alternatively, mixed description QM/MM uses quantum mechanics for the titratable solute molecule and a cluster of solvent molecules in the direct neighborhood and classical molecular mechanics for the more distant solvent molecules.52−57 However, a large number of these works consider titratable residues in watersoluble proteins. Received: May 1, 2016 Published: June 16, 2016 3360

DOI: 10.1021/acs.jctc.6b00446 J. Chem. Theory Comput. 2016, 12, 3360−3369

Article

Journal of Chemical Theory and Computation Table 1. Proton Solvation Energies (kcal mol−1)a

Many approaches that use an implicit representation of the solvent molecules for pKa computations in water can still be considered to be ab initio methods as long as they use no unphysical adjustable parameters such as a linear regression scheme. The action of solvent molecules on solute molecules can for instance be described by Langevin dipoles.58−61 However, the majority of implicit solvent computations use a dielectric continuum to compute solvation energies of small solute molecules in water, a method pioneered by Orttung62 and considerably refined by Honig and co-workers63,64 solving the Poisson equation. An interesting alternative to describe the dielectric continuum is the generalized Born approximation and variants of this approach65−69 instead of solving the Poisson equation. The majority of ab initio pKa computations model the solvent by the dielectric continuum solving the Poisson equation, whereas the energetics of the titratable molecule is evaluated by QC. The first computations of this type were performed by Jorgensen and co-workers70,71 followed by many more contributions over the years. 33,72−87 In several approaches, the solute is dressed by a few explicit water molecules to model the environment of the solute molecule more faithfully.88−90 This involves however some arbitrariness in the number of water molecules and the adopted solute−water conformations. Occasionally, the value of the proton solvation energy in water was avoided by referring to a specific titratable molecule with a known pKa value.91−95 Another issue is pKa values of molecules in different solvents such as acetonitrile (MeCN), dimethyl-sulfoxide (DMSO), and methanol (MeOH). Water is the protic solvent per se because all hydrogens are polar. MeOH involves polar and nonpolar hydrogens and is therefore a mixed solvent (neither protic nor aprotic), whereas MECN and DMSO possess only nonpolar hydrogens and are therefore aprotic solvents. Except for DMSO,96−106 there are fewer molecules with experimentally known pKa values for MeOH98,107,108 and MeCN,109−114 and there is less detailed theoretical work to compute or predict pKa values for those solvents than for water.27,33,68,85,86,91,94,95,115,116 Hence, there is need to close this gap with a cost-efficient ab initio method. We suggest here a method that transforms the pKa value of a titratable residue from one solvent to another by computing the ESE of the protonated and deprotonated molecular species. The ESEs of small molecules and ions in water have been computed reliably for many years.12,15,16,61,64 However, to compute pKa values by ab initio methods, the proton solvation energy of the corresponding solvent must be available. For water, the proton solvation free energy117−119 of −265.9 kcal mol−1 was only recently recognized as the consensus value.19 For other solvents, the precise value of the proton solvation energy is still under debate. An overview of measured and computed proton solvation free energies is given in Table 1. To transform pKa values computationally between different solvents, a crucial point is the knowledge of the proton solvation energy in these solvents. We recently determined the proton solvation energies in MeCN, DMSO, and MeOH by computing the pKa values of a selection of titratable molecules with experimentally known pKa values in at least two of the three considered solvents86,87 (Table 1). The proton solvation energy was determined as a parameter to match the computed pKa values with the measured pKa values. The same procedure was also applied to determine the proton solvation energy in water, yielding

MeCN −255.2b −260.2d −254.2e −254.7f −255.186 −252.3120 −253.2124

DMSO measured −270.5,b −268.6c −273.3d −268.55e computed −266.0f −266.386 −273.2,121 −267.0122 −267.6,120 −261.1124

MeOH −263.8b −263.5d −261.9e −265.4f −265.586 −263.4123 −253.6125

a

The table is updated from Table I in refs 86 and 87 including corrected references. Present results are in bold digits. bBased on transfer free energies from water (using 265.9 kcal mol−1) to other solvents using the TATB approach.126 In earlier work from the same laboratory, the proton solvation free energies were determined by the same method, yielding 254.8, 270.5, and 263.4 kcal mol−1 for MeCN, DMSO, and MeOH, respectively.127 cBased on solvation energies of ions in water and DMSO.128 dBased on the CSB approximation.116 e Based on absolute electrode potentials.129 fBased on the present study using for the proton solvation energy in water −266.3 kcal mol−1 and ignoring the outlier pKa value of 2-chloroaniline in water and MeOH (see Table 3).

−266.3 kcal mol−1, which is very close to the consensus value of 265.9 kcal mol−1, demonstrating the appropriateness of the computational procedure. The error margin86,87 of the proton solvation energies is, for all three solvents, at or below 1 kcal mol−1 corresponding to 0.735 pH units or less at room temperature. The effectiveness of the electrostatic transform procedure in converting pKa values from one solvent to another is demonstrated for 30 compounds belonging to 10 different molecular families, as depicted in Figure 1. They cover a wide range of small organic molecules with variable protonation states. Four different solvents, water, MeCN, DMSO, and MeOH, are considered in this study. However, the electrostatic transform procedure can be applied to any other solvent where the proton solvation energy is known. The known reference

Figure 1. Ten different molecular families involving a total of 30 different titratable molecules (Table 3) displayed in the protonated state were considered for pKa computations using the electrostatic transform method. The number of titratable molecules of a specific family is given in brackets behind the one-letter family name. Titratable oxygens, nitrogens, and sulfurs are highlighted in red, blue, and yellow, respectively. 3361


Article

Journal of Chemical Theory and Computation

solvation energy cancel in the difference of solvation energies of protonated and deprotonated species. Hence, it is sufficient to evaluate only the electrostatic contributions to the solvation energies. All computations are performed with geometryoptimized molecular structures by energy minimization using the density functional theory (DFT) program Jaguar version 8130 with the B3LYP functional, double-ζ with basis set 631G**. Evaluating the ESE of the molecules requires atomic partial charges that are determined by using the restrained electrostatic potential (RESP)131,132 method. The RESP procedure uses the electrostatic potential in the neighborhood of the molecule, which is generated by the nuclear charges and the electronic wave function from a QC computation in the gas phase with the B3LYP functional, double-ζ with basis set 6-31G*. This procedure assumes that gas-phase conformation and atomic partial charges are appropriate. The method works for relatively rigid or small molecules but can be generalized also to larger flexible molecules if one considers appropriate rigid substructures of such molecules. We have generated the atomic partial charges under vacuum conditions. This is done to avoid the electron leakage effect.133,134 The electronic wave function of a molecule embedded in a dielectric medium representing the solvent expands in this medium, an effect called electron leakage. This effect is particularly large in the deprotonated state with an excess negative charge. In reality, there is no dielectric medium, but the electron density of solvent molecules, which repels the electrons of the solute molecule, results in the opposite effect. Even in vacuum, the electronic wave function of a molecule can expand more when embedded in a solvent. To minimize this effect, the atomic partial charges are computed using a low-level electronic wave function: double-ζ with basis set 6-31G**. The ESE terms Gs(BH, j) and Gs(B, j) are evaluated by solving the Poisson equation with the program SOLVATE, from the program suit MEAD.135,136 The investigated solvents are MeOH, DMSO, MeCN, and water, which are described using the dielectric continuum with ε = 32.75, 46.7, 37.5, and 80, respectively, and, inside the volume of the solute, ε = 1. To account for the proper asymptotic electrostatic potential, a twostep focusing scheme is used to solve the Poisson equation with cubic lattices of 1893 grid points with the molecule in the center. The resolutions of the first coarse-grained and the second fine-grained grids are 0.4 and 0.1 Å, respectively. The fine-grained grid is still large enough to contain the complete molecule. This scheme has been tested to yield ESE differences of protonated and deprotonated molecular species with an accuracy of 0.1 kcal mol−1 corresponding to less than 0.1 pH units. The boundary between the solute and the dielectric medium of the solvent is defined by the so-called solvent-accessible surface area (SASA). The SASA is generated by rolling a spherical probe with the solvent-dependent radius over the solute molecule assembled by the joint volumes of atomic spheres. The probe sphere radii are 1.4, 2.05, 2.23, and 2.41 Å for water, MeOH, DMSO, and MeCN, respectively. Because hydrogen atoms have the smallest atomic radius, they can come closest to other not covalently linked atoms. This is particularly important for solvent−solute interactions if the solvent molecule possesses polar hydrogens. To account for such interactions in continuum electrostatic approaches, where the solvent molecule is not explicitly modeled, the solute atomic radii are adjusted and vary with the solvent type.137 They are

pKa value, which is used as the starting point of the electrostatic transform procedure, can be a measured, computed, or predicted by an empirical method. However, the precision of the pKa value that is transformed from the pKa value in another solvent depends on the reliability of this reference pKa value. The method can also be used to check the consistency of pKa values that a molecule has in different solvents. The advantage of the method is that no high-level electronic energy and no costly evaluation of vibrational energies need to be evaluated by QC methods, as is necessary for ab initio computations of pKa values.

2. METHOD 2.1. Theoretical Framework. The key equation to relate the pKa values of a molecule in two different solvents (i and j) that can be in the protonated (BH) or the deprotonated (B) form is [pK a(i) − pK a(j)]kBT ln(10) = ΔGsolv (i) + G H(i) − ΔGsolv (j) − G H(j)

(1)

The terms on the right side of eq 1 are the difference of solvation energies between the BH and B molecular species ΔGsolv (j) = Gsolv (BH, j) − Gsolv (B−, j)

(2)

and the proton solvation energy GH(j) in solvent j. The factor kBT ln(10) converts dimensionless pKa values to energies. At room temperature (T = 298 K), the conversion factor is kBT ln(10) = 1.363 kcal mol−1. Using eqs 1 and 2, we can transform the experimental pKa value, pKexp a,k (j), of a specific molecule k in a solvent j to the pKa value, pKth a,k(j → i), of the same molecule in second solvent i according to pK a,thk(j → i) −1 solv = pK a,exp k (j) + [kBT ln(10)] [ΔGk (i) + G H(i)

− ΔGksolv (j) − G H(j)]

(3)

However, the known pKa value, pKexp a,k (j), can also be a theoretical estimate or a predicted value from other sources. Note that this set of equations is valid only if the binding energy of a proton to a solute molecule does not depend on the environment of the molecule. The success of the method demonstrates that this is indeed the case. In fact, the environment can even be absent (i.e., applying vacuum conditions), and the proton binding energy to the solute is virtually the same as in a solvent. This has been demonstrated before using ab initio methods to compute pKa values,79,86 where the covalent interaction of the proton with the solute was computed in vacuum. The van der Waals interactions between solvent molecules and protonated and deprotonated species of a specific solute molecule are very similar because the two molecular species differ only by one proton. Thus, the only interactions that remain in the difference of pKa values for different solvents are the electrostatic interactions and the proton solvation energies. This is why in the following only the electrostatic solvation energies ΔGsolv k (i) and ΔGsolv k (j) of a titratable molecule k in the solvents i and j are evaluated in eq 3. The electrostatic solvation energies are of course only part of the absolute solvation energies, which are however not needed for pKa computations. 2.2. Electrostatic Solvation Free Energy Computation. It is assumed that all nonelectrostatic contributions to the 3362


Article


3. RESULTS AND DISCUSSION In this work, an electrostatic transform method is introduced to compute macroscopic pKa values of organic molecules in water, MeCN, DMSO, and MeOH from known pKa values of the same molecule in another solvent. To demonstrate the precision and the generality of use, 30 different compounds involving titratable oxygen, nitrogen, or sulfur atoms are considered where the measured pKa values are known for at least two different solvents (Figure 1). In total, 77 measured pKa values in different solvents are considered for which 138 different pKa values are computed by the electrostatic transform method. The same partial atomic charges, computed in vacuum, are used for all of the investigated solvents. The root-mean-square deviations between calculated and available measured pKa values (pKa-RMSDs) are, without the outlier, equal to 0.64, 0.62, 0.66, and 0.58 for water, MeCN, DMSO, and MeOH, considering 29, 17, 18, and 11 different titratable compounds, respectively (see Table 3 and Figure 2). However, considering also the outlier (2-chloroaniline, family I, no. 23), the pKa-

optimized as described in ref 86 and are listed in Table 2. The smallest solute atomic radii are used for water, which is the Table 2. Optimized Atomic Radii for MeCN, MeOH, Water, and DMSO Solvents According to Reference 86 atomic radii (Å) atom

MeCN

MeOH

water

DMSO

H C O/N S

1.00 1.92 1.79

1.00 1.68 1.57 (2.00)

1.00 1.50 1.40 2.00

1.00 1.87 1.75 2.18

protic solvent per se. The solute radii are larger for MeOH, which is neither protic nor aprotic, and largest for DMSO and MeCN, which are purely aprotic solvents possessing no polar hydrogen atom. The proton solvation energies G(H+, j) in a specific solvent j entering the pKa transformation relation, eq 1, are −266.3, −255.1, −266.4, and −265.5 kcal mol−1 for water, MeCN, DMSO, and MeOH, respectively.86,87

Table 3. Measured (exp) and Computed (theo) pKa Values and Their Differences (diff = exp − theo) of 30 Considered Organic Compounds (see Figure 1) in Water, MeCN, DMSO, and MeOHa solvents: no 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

compounds A A A A B C C D D E E E E F F G G H H H H H I I I I J J J J

3,4-methoxybenzoic acid adipic acid benzoic acid fumaric acid phenol benzohydroxamic acid methylhydroxamic acid 3-methylbenzenethiol thiophenol ammonia butylamine ethylamine methylamine diethylamine dimethylamine triethylamine trimethylamine 4-pyrrolidinopyridine imidazole piperidine pyridine pyrrolidine 2-chloroaniline 3-chloroaniline 4-methylaniline aniline methylsulfonamide phenyl sulfonamide Piloty’s acid saccharine pKa-RMSDb

water exp

theo

4.44 4.83 4.41 4.22 4.22 4.44 3.03 3.56 10.00 10.38 8.80 8.30 8.50 9.92 6.66 6.62 6.62 6.37 9.25 10.10 10.80 11.36 10.70 10.73 10.62 11.48 10.98 11.64 10.73 11.48 10.78 11.16 9.80 10.13 9.70 9.11 7.05 6.90 11.12 11.44 5.25 4.73 11.31 12.24 2.66 4.52 3.70 4.05 4.58 5.02 4.60 3.82 10.87 11.20 10.10 10.29 9.29 10.74 1.80 0.49 0.72 (0.64)

acetonitrile diff −0.39 0.19 −0.22 −0.53 −0.38 0.50 −1.22 0.04 0.25 −0.85 −0.56 −0.03 −0.86 −0.66 −0.75 −0.38 −0.33 0.59 −0.15 −0.32 0.52 −0.93 −1.86 −0.35 −0.44 0.78 −0.33 −0.19 1.45 1.31

dimethyl sulfoxide

exp

theo

diff

20.30 20.10 18.60 26.60

20.49 20.20 18.51 26.19

−0.19 −0.10 0.09 0.41

16.50

17.26

−0.76

18.40 18.37 18.75 18.73 18.46 17.61 18.20

18.67 17.54 18.09 18.50 18.08 17.28 18.79

−0.27 0.83 0.66 0.23 0.38 0.33 −0.59

18.92 12.33 19.58

18.90 12.14 18.49

0.02 0.19 1.09

10.56

11.48

−0.92

14.60 0.62

15.91

−1.31

exp

theo

dev

11.40

11.01

0.39

11.10

11.17

−0.17

18.00 13.70 16.00 10.55 10.28 10.50

16.93 14.20 14.78 10.59 10.18 9.83

1.07 −0.50 1.22 −0.04 −0.10 0.67

6.40 10.85 3.50 11.06

6.55 10.60 4.80 11.21

0.15 0.25 −1.30 −0.15

2.98 4.48 3.82 17.50 16.10 15.40

2.63 4.04 3.68 17.17 15.91 13.95

0.35 0.44 0.14 0.33 0.19 −1.45

0.66

methanol exp

theo

diff

9.40 8.00 14.20

8.91 7.56 15.30

0.49 0.44 −1.10

8.95 10.78 11.47 11.00 11.00

9.31 9.84 10.91 10.70 10.97

−0.36 0.94 0.56 0.30 0.03

11.20

10.68

0.52

11.07 5.40

11.02 4.81

0.05 0.59

3.71

1.85

1.86

0.77 (0.58)

a The total number of measured pKa values listed in the table is 77. The number of pKa values computed by the electrostatic transform method is 138. The computed pKa values given in the table for one solvent are the average of the results of the electrostatic transform method from all other solvents with known experimental pKa values. Computed pKa values with a deviation larger than 1.5 pH units are considered to be outliers and marked in bold digits. The only outlier is 2-chloroaniline, family I, no. 23. bIn parentheses, the pKa-RMSD values are given without outlier no. 23: 2chloroaniline.

3363


Article


Figure 2. Correlation diagram of measured and computed pKa values of 30 considered compounds (Figure 1 and Table 3) in water (A), MeCN (B), DMSO (C), and MeOH (D). The number of compounds where the measured pKa values are available is given in brackets behind the pKa-RMSD. The computed pKa values are averages of the pKa values obtained by the electrostatic transform method from the measured pKa values in other solvents. Deviations of more than 1.5 units between the computed and measured pKa values are considered to be outliers and labeled by an open circle (○) in (A) and (D). The pKa-RMSD values without the outlier (2-chloroaniline, family I, no. 23) are given in brackets.

RMSDs increase only moderately to 0.72 and 0.77 in water and MeOH, respectively. There is no outlier for MeCN and DMSO. Generally, there are two main reasons that outliers may occur using the pKa transform method: (1) there are specific solute−solvent conformations that are not considered appropriately by continuum electrostatics and (2) there is a systematic error in the measured pKa value. A third reason could be that the atomic partial charges computed for vacuum conditions are not appropriate in the solvent. The latter should most probably occur for water because it is the most polar one of the considered solvents. But, because the pKa-RMSD value for water is practically as low as that for the other three solvents, this error source cannot be significant. Because the pKa-RMSDs obtained with the electrostatic transform method are in the same range as the values obtained with a complete ab initio method combining QC with electrostatics,79,86,87 the QC part in these approaches is likely more reliable than the electrostatic part in such computations. Deviations between measured and calculated pKa values are larger than 1.5 pH units only for 2-chloroaniline (compound 23). For this compound, the measured pKa values are only available for water and methanol. Hence, upon transforming the pKa value from one solvent to the other, the deviation is of the same absolute value (1.86 pH units) but of different sign (Table 3). Generally, the sum of deviations between computed and measured pKa values of a specific compound in different solvents, listed in Table 3, vanishes. However, because of round-off errors, a value of ±0.01 is possible. The proton solvation energies GH(j) in solvent j have been carefully evaluated for water, MeCN, DMSO, and MeOH by matching the computed and measured pKa values for molecules, which were selected to yield most reliable computed pKa values. However, it is interesting to find out how much these proton solvation energies may differ if we use the present much larger data set of measured and computed pKa values, although the molecules were selected to cover the large space of different titratable molecules in a representative way (Supporting Information). For this purpose, we minimize a quadratic form in the solvation energies with variable proton solvation energies GH(j) for a specific solvent j. For the sake of simplicity, the experimental pKa values are transformed to energies ΔGkexp(i)

= kBT ln(10)

pK a,exp k ( i)

Note that approximately exp solv bound ΔGk (i) = ΔGk (i) + G H(i) − ΔGkbound is ,H , where ΔGk,H the gas-phase free energy difference between the protonated molecule and the deprotonated molecule with a free proton in the gas phase. The resulting quadratic form is a function of the four proton solvation energies GH(j), where j = 1, ..., 4, and reads 30

1 F(G H(1)... G H(4)) = ∑ 4 k=1

4

∑

nk , i

i,j=1 i≠j

nk , j(ΔΔGk (i , j) − G H(i) + G H(j))2

(5)

where ΔΔGk (i , j) = ΔGkexp(i) − ΔGksolv (i) − [ΔGkexp(j) − ΔGksolv (j)] = −ΔΔGk (j , i)

(6)

is approximately equal to GH(i) − GH(j). The factors nk,j consider that experimental pKa values are not available for all combinations of the 30 molecules and the 4 solvents. They are unity if the experimental pKa value of molecule k in the solvent j is available and 0 if not. The value of the quadratic form does not change if a constant energy is added to each proton solvation energy. Hence, upon minimizing the quadratic form by varying GH(i), we can obtain only relative energy values. To obtain absolute values, we fix the proton solvation energy of water GH(w) = GH(4), whose value is the supposed to be the most precise. Taking the partial derivatives of the function F, eq 5, with respect to the other three proton solvation energies, GH(m), where m = 1, ..., 3, and setting them to 0 yields three linear equations in the unknown GH(m) 3

∑ [Ni ,mG H(i)] − NmG H(m) = ΔΔG(m) − N4,mG H(4) i=1 i≠m

m = 1, ..., 3

(4)

(7)

where 3364


Article

Journal of Chemical Theory and Computation 30

ΔΔG(m) =

4

molecule. Specific solvent-dependent atomic radii are used to define the effective volume of the solute molecules to account for the amount of polar hydrogen atoms of the solvent molecules. Because hydrogen atoms of the solvent can come closer to atoms of the solute molecule than other atoms of the solvent, the solute atomic radii diminish proportional to the fraction of polar hydrogens in solvent molecules. Finally, the socalled SASA separating the dielectric continuum of the solvent from the solute is obtained by rolling with a solvent sphere over the solute volume using the rolling sphere surface patches that are closest to the solute. Another very essential quantity is the proton solvation energy, which was evaluated in preceding work for water, MeCN, DMSO, and MeOH, matching measured and ab initio computed pKa values for a specific selection of molecules. We note that in contrast to other work no explicit solvent molecules needed to be considered in these computations. The pKa-RMSD of the computed and measured pKa values is about 0.7 pH units or even below if the few outliers are omitted. This precision is close to alternative computationally more expensive ab initio methods and can only be outperformed by empirical methods that involve a large number of molecular-family-dependent adjustable parameters. We consider 0.7 pH units as the limiting precision for computational procedures that use a combination of QC and dielectric continuum theory without empirical fetch factors. For higher precision, the explicit solvent structure in the neighborhood of the solute molecule as well as the nonelectrostatic solvation energies will be needed. Furthermore, a separate evaluation of the quantum chemical and electrostatic energies used in simplified ab initio methods for pK a computations79,86,87 is an approximation that may be too crude to yield a precision better than 0.7 pH units. Methods that avoid such a separate evaluation are under development. The computed pKa values can be used to determine the relative proton solvation energies between the considered solvents. Absolute values are obtained if the absolute value of the solvation energy is given for one solvent. Choosing water as a reference, we obtained proton solvation energies for MeCN, DMSO, and MeOH, which are very close to the values we have obtained recently by an ab initio method using QC and electrostatic energy computations.86,87 This demonstrates the reliability of the proton solvation energies and at the same time also the reliability of the electrostatic transform method.

∑ ∑ nk ,i nk ,mΔΔGk(i , m) k=1 i=1 i≠m

(8)

and 30

Ni , m =

∑ nk ,i nk ,m k=1

(9)

is the number of molecules for which in the two solvents i and m experimental pKa values are available and 4

Nm =

∑ Ni ,m i=1 i≠m

(10)

is the number of known pKa values for each considered molecule, where the pKa value in solvent m is also known. The results of the computations of proton solvation energies are −254.68 (−254.71), −266.01 (−266.03), and −265.44 (−265.52) kcal mol−1 for MeCN, DMSO, and MeOH, respectively. In brackets are the values obtained including the outlier pKa value of 2-chloroaniline, which are nearly the same. These proton solvation energies were obtained using 266.3 kcal mol−1 for the proton solvation in water. Alternatively, one could use the consensus value of proton solvation in water as reference, which is 265.9 kcal mol−1, yielding the values −254.28, −265.61, and −265.04 kcal mol−1 for MeCN, DMSO, and MeOH, respectively. These values are very close to the proton solvation energies, which were recently computed by an ab initio method using a combination of QC and continuum electrostatics.86,87 In this approach, we considered molecules that were selected to yield the most reliable values of computed pKa’s. The largest difference of 0.5 kcal mol−1 occurs for MeCN. For DMSO and MeOH, the differences are only 0.2 and 0.0 kcal mol−1. Re-evaluating the pKa-RMSDs with the new proton solvation energies from Table 1, we obtain 0.62, 0.69, 0.58, and 0.47 for water, MeCN, DMSO, and MeOH, respectively. For MeCN, the pKa-RMSD is slightly larger, whereas it is slightly smaller for the other three solvents (see the last row in Table 3) than using the proton solvation energies from refs 86 and 87.

4. SUMMARY We have introduced an electrostatic transform method to compute pKa values for a large variety of titratable compounds in different solvents, which can be applied if the pKa value in one solvent is known either by experiment or by computation. The computation demands only a moderate amount of CPU time because no vibrational energy and electronic ground-state energy at a high level of precision need to be computed. An equilibrium geometry of the molecule is needed, usually obtained by QC geometry optimization in vacuum with double-ζ basis set. For this geometry, atomic partial charges are computed by matching the electrostatic potential generated by these charges with the electrostatic potential from the nuclear charges and electronic wave function of the molecule. Only electrostatic contributions of solvation energy matter for the pKa computations. The nonelectrostatic solvation energies essentially cancel in the difference between protonated and deprotonated species of the same molecule. Very essential for the accuracy of the computations is a proper definition of the molecular surface of the solute used for the evaluation of the ESE of the protonated and deprotonated species of the

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.6b00446. A list of all measured and computed pKa (Tables S1−S4) values as well as the atomic coordinates and charges of all 30 considered molecules (Table S5) (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors are grateful for useful discussions with Dr Art Bochevarov. This work was financially supported by the 3365


Article


(17) Lee, A. C.; Crippen, G. M. Predicting pKa. J. Chem. Inf. Model. 2009, 49, 2013−2033. (18) Ho, J.; Coote, M. L. A universal approach for continuum solvent pKa calculations are we there yet? Theor. Chem. Acc. 2010, 125, 3−21. (19) Alongi, K. S.; Shields, G. C. Theoretical calculations of acid dissociation constants: a review article. In Annual Reports in Computational Chemistry; Wheeler, R. A., Ed.; Elsevier: Amsterdam, 2010; Vol. 6, pp 113−138. (20) Ho, J.; Coote, M. L. First-principles prediction of acidities in the gas and solution phase. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2011, 1, 649−660. (21) Meloun, M.; Bordovská, S. Benchmarking and validating algorithms that estimate pKa values of drugs based on their molecular structures. Anal. Bioanal. Chem. 2007, 389, 1267−1281. (22) Rupp, M.; Köerner, R.; Tetko, I. V. Estimation of acid dissociation constants using graph kernels. Mol. Inf. 2010, 29, 731− 741. (23) Ščavničar, A.; Balaban, A. T.; Pompe, M. Application of variable anti-connectivity index to active sites. Modelling pKa values of aliphatic monocarboxylic acids. SAR QSAR Environ. Res. 2013, 24, 553−563. (24) Tehan, B. G.; Lloyd, E. J.; Wong, M. G.; Pitt, W. R.; Montana, J. G.; Manallack, D. J.; Gancia, E. Estimation of pKa using semiempirical molecular orbital methods. Part 1: application to phenols and carboxylic acids. Quant. Struct.-Act. Rel. 2002, 21, 457−472. (25) Tehan, B. G.; Lloyd, E. J.; Wong, M. G.; Pitt, W. R.; Gancia, E.; Manallack, D. J. Estimation of pKa using semiempirical molecular orbital methods. Part 2: application to amines, anilines and various nitrogen containing heterocyclic compounds. Quant. Struct.-Act. Rel. 2002, 21, 473−485. (26) Krieger, E.; Nielsen, J. E.; Spronk, C. A. E. M.; Vriend, G. Fast empirical pKa prediction by Ewald summation. J. Mol. Graphics Modell. 2006, 25, 481−486. (27) Jover, J.; Bosque, R.; Sales, J. Neural network based QSPR study for predicting pKa of phenols in different solvents. QSAR Comb. Sci. 2007, 26, 385−397. (28) Ugur, I.; Marion, A.; Parant, S.; Jensen, J. H.; Monard, G. Rationalization of the pKa values of alcohols and thiols using atomic charge descriptors and its application to the prediction of amino acid pKa’s. J. Chem. Inf. Model. 2014, 54, 2200−2213. (29) Geidl, S.; Svobodová Vareková, R.; Bendová, V.; Petrusek, L.; Ionescu, C.-M.; Jurka, Z.; Abagyan, R.; Koča, J. J. How does the methodology of 3D structure preparation influence the quality of pKa prediction? Chem. Inf. Model. 2015, 55, 1088−1097. (30) Rupp, M.; Ramakrishnan, R.; von Lilienfeld, O. A. Machine learning for quantum mechanical properties of atoms in molecules. J. Phys. Chem. Lett. 2015, 6, 3309−3313. (31) Parthasarathi, R.; Padmanabhan, J.; Elango, M.; Chitra, K.; Subramanian, V.; Chattaraj, P. K. pKa prediction using group philicity. J. Phys. Chem. A 2006, 110, 6540−6544. (32) Klicić, J. J.; Friesner, R. A.; Liu, S.-Y.; Guida, W. C. Accurate prediction of acidity constants in aqueous solution via density functional theory and self-consistent reaction field methods. J. Phys. Chem. A 2002, 106, 1327−1335. (33) Chipman, D. M. Computation of pKa from dielectric continuum theory. J. Phys. Chem. A 2002, 106, 7413−7422. (34) Klamt, A.; Eckert, F.; Diedenhofen, M.; Beck, M. E. First principles calculations of aqueous pKa values for organic and inorganic acids using COSMO−RS reveal an inconsistency in the slope of the pKa scale. J. Phys. Chem. A 2003, 107, 9380−9386. (35) Eckert, F.; Klamt, A. Accurate prediction of basicity in aqueous solution with COSMO-RS. J. Comput. Chem. 2006, 27, 11−19. (36) Lu, H.; Chen, X.; Zhan, C.-G. First-principles calculation of pKa for cocaine, nicotine, neurotransmitters, and anilines in aqueous solution. J. Phys. Chem. B 2007, 111, 10599−10605. (37) Zhang, S.; Baker, J.; Pulay, P. A reliable and efficient first principles-based method for predicting pKa values. 1. Methodology. J. Phys. Chem. A 2010, 114, 425−431.

German Science Foundation (DFG) Project C2 in the Collaborative Research Center (Sfb1078).

■

ABBREVIATIONS DFT, density functional theory; DMSO, dimethyl sulfoxide; ESE, electrostatic solvation energy; MD, molecular dynamics; MeCN, acetonitrile; MeOH, methanol; pKa-RMSD, root-meansquare deviation of pKa values; QC, quantum chemistry; RESP, restrained electrostatic potential; SASA, solvent-accessible surface area

■

REFERENCES

(1) Katritzky, A. R.; Kuanar, M.; Slavov, S.; Hall, C. D.; Karelson, M.; Kahn, I.; Dobchev, D. A. Quantitative correlation of physical and chemical properties with chemical structure: utility for prediction. Chem. Rev. (Washington, DC, U.S.) 2010, 110, 5714−5789. (2) Mamy, L.; Patureau, D.; Barriuso, E.; Bedos, C.; Bessac, F.; Louchart, X.; Martin-Laurent, F.; Miege, C.; Benoit, P. Prediction of the fate of organic compounds in the environment from their molecular properties: a review. Crit. Rev. Environ. Sci. Technol. 2015, 45, 1277−1377. (3) Nieto-Draghi, C.; Fayet, G.; Creton, B.; Rozanska, X.; Rotureau, P.; de Hemptinne, J.-C.; Ungerer, P.; Rousseau, B.; Adamo, C. A general guidebook for the theoretical prediction of physicochemical properties of chemicals for regulatory purposes. Chem. Rev. (Washington, DC, U.S.) 2015, 115, 13093−13164. (4) Manallack, D. T. The pKa distribution of drugs application to drug discovery. Perspect. Med. Chem. 2007, 1, 25−38. (5) Kerns, E. H.; Di, L. Drug-like Properties: Concepts, Structure Design and Methods: From ADME to Toxicity Optimization; Academic Press: London, 2008. (6) Manallack, D. T.; Prankerd, R. J.; Yuriev, E.; Oprea, T. I.; Chalmers, D. K. The significance of acid/base properties in drug discovery. Chem. Soc. Rev. 2013, 42, 485−496. (7) Diaz, D.; Ford, K. A.; Hartley, D. P.; Harstad, E. B.; Cain, G. R.; Achilles-Poon, K.; Nguyen, T.; Peng, J.; Zheng, Z.; Merchant, M.; Sutherlin, D. P.; Gaudino, J. J.; Kaus, R.; Lewin-Koh, S. C.; Choo, E. F.; Liederer, B. M.; Dambach, D. M. Pharmacokinetic drivers of toxicity for basic molecules: Strategy to lower pKa results in decreased tissue exposure and toxicity for a small molecule Met inhibitor. Toxicol. Appl. Pharmacol. 2013, 266, 86−94. (8) Charifson, P. S.; Walters, W. P. Acidic and basic drugs in medicinal chemistry: a perspective. J. Med. Chem. 2014, 57, 9701− 9717. (9) Jinich, A.; Rappoport, D.; Dunn, I.; Sanchez-Lengeling, B.; Olivares-Amaya, R.; Noor, E.; Bar Even, A.; Aspuru-Guzik, A. Quantum chemical approach to estimating the thermodynamics of metabolic reactions. Sci. Rep. 2014, 4, No. 7022. (10) Kollman, P. Free energy calculations: applications to chemical and biochemical phenomena. Chem. Rev. (Washington, DC, U.S.) 1993, 93, 2395−2417. (11) Tomasi, J.; Persico, M. Molecular interactions in solution: an overview of methods based on continuous distributions of the solvent. Chem. Rev. (Washington, DC, U.S.) 1994, 94, 2027−2094. (12) Cramer, C. J.; Truhlar, D. G. Implicit solvation models: equilibria, structure, spectra, and dynamics. Chem. Rev. (Washington, DC, U.S.) 1999, 99, 2161−2200. (13) Orozco, M.; Luque, F. J. Theoretical methods for the description of the solvent effect in biomolecular systems. Chem. Rev. (Washington, DC, U.S.) 2000, 100, 4187−4226. (14) Simonson, T. Electrostatics and dynamics of proteins. Rep. Prog. Phys. 2003, 66, 737−787. (15) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum mechanical continuum solvation models. Chem. Rev. (Washington, DC, U.S.) 2005, 105, 2999−3094. (16) Cramer, C. J.; Truhlar, D. G. A universal approach to solvation modeling. Acc. Chem. Res. 2008, 41, 760−768. 3366


Article

Journal of Chemical Theory and Computation (38) Zhang, S. A reliable and efficient first principles-based method for predicting pKa values. 4. Organic bases. J. Comput. Chem. 2012, 33, 2469−2482. (39) Matsui, T.; Baba, T.; Kamiya, K.; Shigeta, Y. An accurate density functional theory based estimation of pKa values of polar residues combined with experimental data: from amino acids to minimal proteins. Phys. Chem. Chem. Phys. 2012, 14, 4181−4187. (40) Warshel, A. Dynamics of reactions in polar solvents. Semiclassical trajectory studies of electron-transfer and proton-transfer reactions. J. Phys. Chem. 1982, 86, 2218−2224. (41) Warshel, A.; Sussman, F.; King, G. Free energy of charges in solvated proteins: microscopic calculations using a reversible charging process. Biochemistry 1986, 25, 8368−8372. (42) Bürgi, R.; Kollman, P. A.; Van Gunsteren, W. F. Simulating proteins at constant pH: an approach combining molecular dynamics and Monte Carlo simulation. Proteins: Struct., Funct., Bioinf. 2002, 47, 469−480. (43) Baptista, A. M.; Teixeira, V. H.; Soares, C. M. Constant-pH molecular dynamics using stochastic titration. J. Chem. Phys. 2002, 117, 4184−4200. (44) Simonson, T.; Carlsson, J.; Case, D. A. Proton binding to proteins: pKa calculations with explicit and implicit solvent models. J. Am. Chem. Soc. 2004, 126, 4167−4180. (45) Arthur, E. J.; Yesselman, J. D.; Brooks, C. L. Predicting extreme pKa shifts in staphylococcal nuclease mutants with constant pH molecular dynamics. Proteins: Struct., Funct., Bioinf. 2011, 79, 3276− 3286. (46) Bonthuis, D. J.; Netz, R. R. Beyond the continuum: how molecular solvent structure affects electrostatics and hydrodynamics at solid−electrolyte interfaces. J. Phys. Chem. B 2013, 117, 11397−11413. (47) Lin, Y.-L.; Aleksandrov, A.; Simonson, T.; Roux, B. An overview of electrostatic free energy computations for solutions and proteins. J. Chem. Theory Comput. 2014, 10, 2690−2709. (48) Ivanov, I.; Chen, B.; Raugei, S.; Klein, M. L. Relative pKa values from first-principles molecular dynamics: the case of histidine deprotonation. J. Phys. Chem. B 2006, 110, 6365−6371. (49) Cheng, J.; Liu, X.; VandeVondele, J.; Sulpizi, M.; Sprik, M. Redox potentials and acidity constants from density functional theory based molecular dynamics. Acc. Chem. Res. 2014, 47, 3522−3529. (50) Tummanapelli, A. K.; Vasudevan, S. Estimating successive pKa values of polyprotic acids from ab initio molecular dynamics using metadynamics: the dissociation of phthalic acid and its isomers. Phys. Chem. Chem. Phys. 2015, 17, 6383−6388. (51) Car, R.; Parrinello, M. Unified approach for molecular dynamics and density-functional theory. Phys. Rev. Lett. 1985, 55, 2471−2474. (52) Davies, J. E.; Doltsinis, N. L.; Kirby, A. J.; Roussev, C. D.; Sprik, M. Estimating pKa values for pentaoxyphosphoranes. J. Am. Chem. Soc. 2002, 124, 6594−6599. (53) Riccardi, D.; Schaefer, P.; Cui, Q. pKa calculations in solution and proteins with QM/MM free energy perturbation simulations: a quantitative test of QM/MM protocols. J. Phys. Chem. B 2005, 109, 17715−17733. (54) Sulpizi, M.; Sprik, M. Acidity constants from vertical energy gaps: density functional theory based molecular dynamics implementation. Phys. Chem. Chem. Phys. 2008, 10, 5238−5249. (55) Ghosh, N.; Cui, Q. pKa of residue 66 in Staphylococal nuclease. I. Insights from QM/MM simulations with conventional sampling. J. Phys. Chem. B 2008, 112, 8387−8397. (56) Kamerlin, S. C. L.; Haranczyk, M.; Warshel, A. Progress in ab initio QM/MM free-energy simulations of electrostatic energies in proteins: accelerated QM/MM studies of pKa, redox reactions and solvation free energies. J. Phys. Chem. B 2009, 113, 1253−1272. (57) Mangold, M.; Rolland, L.; Costanzo, F.; Sprik, M.; Sulpizi, M.; Blumberger, J. Absolute pKa values and solvation structure of amino acids from density functional based molecular dynamics simulation. J. Chem. Theory Comput. 2011, 7, 1951−1961. (58) Russell, S. T.; Warshel, A. Calculations of electrostatic energies in proteins. The energetics of ionized groups in bovine pancreatic trypsin inhibitor. J. Mol. Biol. 1985, 185, 389−404.

(59) Warshel, A.; Aqvist, J. Electrostatic energy and macromolecular function. Annu. Rev. Biophys. Biophys. Chem. 1991, 20, 267−298. (60) Luzhkov, V.; Warshel, A. Microscopic models for quantum mechanical calculations of chemical processes in solutions: LD/ AMPAC and SCAAS/AMPAC calculations of solvation energies. J. Comput. Chem. 1992, 13, 199−213. (61) Florián, J.; Warshel, A. Langevin dipoles model for ab initio calculations of chemical processes in solution: parameterization and application to hydration free energies of neutral and ionic solutes and conformational analysis in aqueous solution. J. Phys. Chem. B 1997, 101, 5583−5595. (62) Orttung, W. H. Direct solution of the Poisson equation for biomolecules of arbitrary shape, polarizability density and charge distribution. Ann. N.Y. Acad. Sci. 1977, 303, 22−37. (63) Gilson, M. K.; Sharp, K. A.; Honig, B. Calculating the electrostatic potential of molecules in solution: method and error assessment. J. Comput. Chem. 1988, 9, 327−335. (64) Sitkoff, D.; Sharp, K. A.; Honig, B. Accurate calculation of hydration free energies using macroscopic solvent models. J. Phys. Chem. 1994, 98, 1978−1988. (65) Jayaram, B.; Liu, Y.; Beveridge, D. L. A modification of the generalized Born theory for improved estimates of solvation energies and pK shifts. J. Chem. Phys. 1998, 109, 1465−1471. (66) Hassan, S. A.; Guarnieri, F.; Mehler, E. L. A general treatment of solvent effects based on screened coulomb potentials. J. Phys. Chem. B 2000, 104, 6478−6489. (67) Sharma, I.; Kaminski, G. A. Calculating pKa values for substituted phenols and hydration energies for other compounds with the first-order fuzzy-border continuum solvation model. J. Comput. Chem. 2012, 33, 2388−2399. (68) Pomogaeva, A.; Chipman, D. M. Hydration energy from a composite method for implicit representation of solvent. J. Chem. Theory Comput. 2014, 10, 211−219. (69) Pomogaeva, A.; Chipman, D. M. Composite method for implicit representation of solvent in dimethyl sulfoxide and acetonitrile. J. Phys. Chem. A 2015, 119, 5173−5180. (70) Jorgensen, W. L.; Briggs, J. M.; Gao, J. A priori calculations of pKa’s for organic compounds in water. The pKa of ethane. J. Am. Chem. Soc. 1987, 109, 6857−6858. (71) Jorgensen, W. L.; Briggs, J. M. A priori pKa calculations and the hydration of organic anions. J. Am. Chem. Soc. 1989, 111, 4190−4197. (72) Lim, C.; Bashford, D.; Karplus, M. Absolute pKa calculations with continuum dielectric methods. J. Phys. Chem. A 1991, 95, 5610− 5620. (73) Schüürmann, G.; Cossi, M.; Barone, V.; Tomasi, J. Prediction of the pKa of carboxylic acids using the ab initio continuum-solvation model PCM-UAHF. J. Phys. Chem. A 1998, 102, 6706−6712. (74) 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. (75) Jang, Y. H.; Sowers, L. C.; Ç aǧin, T.; Goddard, W. A., III First principles calculation of pKa values for 5-substituted uracils. J. Phys. Chem. A 2001, 105, 274−280. (76) Liptak, M. D.; Shields, G. C. Accurate pKa calculations for carboxylic acids using complete basis set and Gaussian-n models combined with CPCM continuum solvation methods. J. Am. Chem. Soc. 2001, 123, 7314−7319. (77) Liptak, M. D.; Gross, K. C.; Seybold, P. J.; Feldgus, S.; Shields, C. G. Absolute pKa determinations for substituted phenols. J. Am. Chem. Soc. 2002, 124, 6421−6427. (78) Saracino, G.; Improta, R.; Barone, V. Absolute pK a determination for carboxylic acids using density functional theory and the polarizable continuum model. Chem. Phys. Lett. 2003, 373, 411−415. (79) Schmidt am Busch, M.; Knapp, E. W. Accurate pK a determination for a heterogeneous set of organic molecules. ChemPhysChem 2004, 5, 1513−1522. 3367


Article

Journal of Chemical Theory and Computation (80) Jensen, J. H.; Li, H.; Robertson, A. D.; Molina, P. A. Prediction and rationalization of protein pKa values using QM and QM/MM methods. J. Phys. Chem. A 2005, 109, 6634−6643. (81) Takano, Y.; Houk, N. Benchmarking the conductor-like polarizable continuum model (CPCM) for aqueous solvation free energies of neutral and ionic organic molecules. J. Chem. Theory Comput. 2005, 1, 70−77. (82) da Silva, G.; Kennedy, E. M.; Dlugogorski, B. Z. Ab initio procedure for aqueous-phase pKa calculation: the acidity of nitrous acid. J. Phys. Chem. A 2006, 110, 11371−11376. (83) Ho, J.; Coote, M. L. pKa calculation of some biologically important carbon acidsan assessment of contemporary theoretical procedures. J. Chem. Theory Comput. 2009, 5, 295−306. (84) Casasnovas, R.; Ortega-Castro, J.; Frau, J.; Donoso, J.; Muñoz, F. Theoretical pKa calculations with continuum model solvents, alternative protocols to thermodynamic cycles. Int. J. Quantum Chem. 2014, 114, 1350−1363. (85) Farrokhpour, H.; Manassir, M. Approach for predicting the standard free energy solvation of H+ and acidity constant in nonaqueous organic solvents. J. Chem. Eng. Data 2014, 59, 3555− 3564. (86) Rossini, E.; Knapp, E. W. Proton solvation in protic and aprotic solvents. J. Comput. Chem. 2016, 37, 1082−1091. (87) Rossini, E.; Knapp, E. W. Erratum: proton solvation in protic and aprotic solvents. J. Comput. Chem. 2016. (88) Pliego, J. R., Jr.; Riveros, J. M. Theoretical calculation of pKa using the cluster−continuum model. J. Phys. Chem. A 2002, 106, 7434−7439. (89) Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. Adding explicit solvent molecules to continuum solvent calculations for the calculation of aqueous acid dissociation constants. J. Phys. Chem. A 2006, 110, 2493−2499. (90) Zhang, S. A reliable and efficient first principles-based method for predicting pKa values. III. Adding explicit water molecules: can the theoretical slope be reproduced and pKa values predicted more accurately? J. Comput. Chem. 2012, 33, 517−526. (91) Almerindo, G. I.; Tondo, D. W.; Pliego, J. R. Ionization of organic acids in dimethyl sulfoxide solution: a theoretical ab initio calculation of the pKa using a new parametrization of the polarizable continuum model. J. Phys. Chem. A 2004, 108, 166−171. (92) Govender, K. K.; Cukrowski, I. Density functional theory and isodesmic reaction based prediction of four stepwise protonation constants, as log KH(n), for nitrilotriacetic acid. The importance of a kind and protonated form of a reference molecule used. J. Phys. Chem. A 2010, 114, 1868−1878. (93) Sastre, S.; Casasnovas, R.; Muñoz, F.; Frau, J. Isodesmic reaction for pKa calculations of common organic molecules. Theor. Chem. Acc. 2013, 132, 1310−1317. (94) Poliak, P. The DFT calculations of pKa values of the cationic acids of aniline and pyridine derivatives in common solvents. Acta Chim. Slovaca 2014, 7, 25−30. (95) Miguel, E. L. M.; Silva, P. L.; Pliego, J. R. Theoretical prediction of pKa in methanol: testing SM8 and SMD models for carboxylic acids, phenols, and amines. J. Phys. Chem. B 2014, 118, 5730−5739. (96) Bordwell, F. G.; Algrim, D. Nitrogen acids. 1. Carboxamides and sulfonamides. J. Org. Chem. 1976, 41, 2507−2508. (97) Olmstead, W. N.; Bordwell, F. G. Ion-pair association constants in dimethyl sulfoxide. J. Org. Chem. 1980, 45, 3299−3305. (98) Bordwell, F. G.; Hughes, D. L. Thiol acidities and thiolate ion reactivities toward butyl chloride in dimethyl sulfoxide solution. The question of curvature in Broensted plots. J. Org. Chem. 1982, 47, 3224−3232. (99) Bordwell, F. G.; McCallum, R. J.; Olmstead, W. N. Acidities and hydrogen bonding of phenols in dimethyl sulfoxide. J. Org. Chem. 1984, 49, 1424−1427. (100) Bordwell, F. G. Equilibrium acidities in dimethyl sulfoxide solution. Acc. Chem. Res. 1988, 21, 456−463. (101) Bordwell, F. G.; Fried, H. E.; Hughes, D. L.; Lynch, T. Y.; Satish, A. V.; Whang, Y. E. Acidities of carboxamides, hydroxamic

acids, carbohydrazides, benzenesulfonamides, and benzenesulfonohydrazides in DMSO solution. J. Org. Chem. 1990, 55, 3330−3336. (102) Bordwell, F. G.; Harrelson, J. A., Jr.; Lynch, T. Y. Homolytic bond dissociation energies for the cleavage of alpha-nitrogen-hydrogen bonds in carboxamides, sulfonamides, and their derivatives. The question of synergism in nitrogen-centered radicals. J. Org. Chem. 1990, 55, 3337−3341. (103) Bordwell, F. G.; Satish, A. V. Acidities of C2 hydrogen atoms in thiazolium cations and reactivities of their conjugate bases. J. Am. Chem. Soc. 1991, 113, 985−990. (104) Crampton, M. R.; Robotham, I. A. Acidities of some substituted ammonium ions in dimethyl sulfoxide. J. Chem. Res., Synop. 1997, 22−23. (105) Ripin, D. H.; Evans, D. A. Chem 206, pKa’s of Inorganic and Oxoacids. http://evans.rc.fas.harvard.edu/pdf/evans_pKa_table.pdf (2005). (106) Jover, J.; Bosque, R.; Sales, J. QSPR prediction of pKa for benzoic acids in different solvents. QSAR Comb. Sci. 2008, 27, 563− 581. (107) Rived, F.; Rosés, M.; Bosch, E. Dissociation constants of neutral and charged acids in methyl alcohol. The acid strength resolution. Anal. Chim. Acta 1998, 374, 309−324. (108) Rived, F.; Canals, I.; Bosch, E.; Rosés, M. Acidity of methanol−water. Anal. Chim. Acta 2001, 439, 315−333. (109) Coetzee, J. F.; Padmanabhan, G. R. Dissociation and homoconjugation of certain phenols in acetonitrile. J. Phys. Chem. 1965, 69, 3193−3196. (110) Kolthoff, I. M.; Chantooni, M. K. Conductometric, potentiometric, and spectrophotometric determination of dissociation constants of substituted benzoic acids in acetonitrile. J. Phys. Chem. 1966, 70, 856−866. (111) Chantooni, M. K.; Kolthoff, I. M. Electronic absorption spectra of ion pairs composed of substituted amine picrates in acetonitrile. J. Am. Chem. Soc. 1968, 90, 3005−3009. (112) Espinosa, S.; Bosch, E.; Rosés, M. Retention of ionizable compounds in high-performance liquid chromatography. 14. Acid−base pK values in acetonitrile−water mobile phases. J. Chromatogr. A 2002, 964, 55−66. (113) Nicoleti, C. R.; Marini, V. G.; Zimmermann, L. M.; Machado, V. G. Anionic chromogenic chemosensors highly selective for fluoride or cyanide based on 4-(4-nitrobenzylideneamine)-phenol. J. Braz. Chem. Soc. 2012, 23, 1488−1500. (114) Raamat, E.; Kaupmees, K.; Ovsjannikov, G.; Trummal, A.; Kütt, A.; Saame, J.; Koppel, I.; Kaljurand, I.; Lipping, L.; Rodima, T.; Pihl, V.; Koppel, I. A.; Leito, I. Acidities of strong neutral Brønsted acids in different media. J. Phys. Org. Chem. 2013, 26, 162−170. (115) Li, J.-N.; Fu, Y.; Liu, L.; Guo, Q.-X. First-principle predictions of basicity of organic amines and phosphines in acetonitrile. Tetrahedron 2006, 62, 11801−11813. (116) Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. Single-ion solvation free energies and the normal hydrogen electrode potential in methanol, acetonitrile, and dimethyl sulfoxide. J. Phys. Chem. B 2007, 111, 408−422. (117) Tissandier, M. D.; Cowen, K. A.; Feng, W. Y.; Gundlach, E.; Cohen, M. H.; Earhart, A. D.; Coe, J. V.; Tuttle, T. R. The proton’s absolute aqueous enthalpy and Gibbs free energy of solvation from cluster-ion solvation data. J. Phys. Chem. A 1998, 102, 7787−7794. (118) Camaioni, D. M.; Schwerdtfeger, C. A. Comment on “accurate experimental values for the free energies of hydration of H+, OH−, and H3O+”. J. Phys. Chem. A 2005, 109, 10795−10797. (119) Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. Aqueous solvation free energies of ions and ion−water clusters based on an accurate value for the absolute aqueous solvation free energy of the proton. J. Phys. Chem. B 2006, 110, 16066−16081. (120) Himmel, D.; Goll, S. K.; Leito, I.; Krossing, I. Anchor points for the unified Brønsted acidity scale: the rCCC model for the calculation of standard Gibbs energies of proton solvation in eleven representative liquid media. Chem.−Eur. J. 2011, 17, 5808−5826. 3368


Article

Journal of Chemical Theory and Computation (121) Westphal, E.; Pliego, J. R. Absolute solvation free energy of Li+ and Na+ ions in dimethyl sulfoxide solution: a theoretical ab initio and cluster-continuum model study. J. Chem. Phys. 2005, 123, No. 074508. (122) Fu, Y.; Liu, L.; Li, R.-Q.; Liu, R.; Guo, Q.-X. First-principle predictions of absolute pKa’s of organic acids in dimethyl sulfoxide solution. J. Am. Chem. Soc. 2004, 126, 814−822. (123) Hwang, S.; Chung, D. S. Calculation of the solvation free energy of the proton in methanol. Bull. Korean Chem. Soc. 2005, 26, 589−593. (124) Carvalho, N. F.; Pliego, J. R., Jr. Cluster-continuum quasichemical theory calculation of the lithium ion solvation in water, acetonitrile and dimethyl sulfoxide: an absolute single-ion solvation free energy scale. Phys. Chem. Chem. Phys. 2015, 17, 26745− 26755. (125) Pliego, J. R.; Miguel, E. L. M. Absolute single-ion solvation free energy scale in methanol determined by the lithium cluster-continuum approach. J. Phys. Chem. B 2013, 117, 5129−5135. (126) Kalidas, C.; Hefter, G.; Marcus, Y. Gibbs energies of transfer of cations from water to mixed aqueous organic solvents. Chem. Rev. (Washington, DC, U.S.) 2000, 100, 819−852. (127) Marcus, Y.; Kamlet, M. J.; Taft, R. W. Linear solvation energy relationships: standard molar Gibbs free energies and enthalpies of transfer of ions from water into nonaqueous solvents. J. Phys. Chem. 1988, 92, 3613−3622. (128) Pliego, J. R., Jr.; Riveros, J. M. Gibbs energy of solvation of organic ions in aqueous and dimethyl sulfoxide solutions. Phys. Chem. Chem. Phys. 2002, 4, 1622−1627. (129) Fawcett, W. R. The ionic work function and its role in estimating absolute electrode potentials. Langmuir 2008, 24, 9868− 9875. (130) Bochevarov, A. D.; Harder, E.; Hughes, T. F.; Greenwood, J. R.; Braden, D. A.; Philipp, D. M.; Rinaldo, D.; Halls, M. D.; Zhang, J.; Friesner, R. A. Jaguar: a high-performance quantum chemistry software program with strengths in life and materials sciences. Int. J. Quantum Chem. 2013, 113, 2110−2142. (131) Besler, B. H.; Merz, K. M.; Kollman, P. A. Atomic charges derived from semiempirical methods. J. Comput. Chem. 1990, 11, 431− 439. (132) Bayly, C. I.; Cieplak, P.; Cornell, W.; Kollman, P. A. A wellbehaved electrostatic potential based method using charge restraints for deriving atomic charges: the RESP model. J. Phys. Chem. 1993, 97, 10269−10280. (133) Richardson, W. H.; Peng, C.; Bashford, D.; Noodleman, L.; Case, D. A. Incorporating solvation effects into density functional theory: calculation of absolute acidities. Int. J. Quantum Chem. 1997, 61, 207−217. (134) Schmidt am Busch, M.; Knapp, E. W. One-electron reduction potential for oxygen- and sulfur-centered organic radicals in protic and aprotic solvents. J. Am. Chem. Soc. 2005, 127, 15730−15737. (135) Bashford, D.; Gerwert, K. Electrostatic calculations of the pKa values of ionizable groups in bacteriorhodopsin. J. Mol. Biol. 1992, 224, 473−486. (136) Bashford, D.; Case, D. A.; Dalvit, C.; Tennant, L.; Wright, P. E. Electrostatic calculations of side-chain pKa values in myoglobin and comparison with NMR data for histidines. Biochemistry 1993, 32, 8045−8056. (137) Marenich, A. V.; Olson, R. M.; Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. Self-consistent reaction field model for aqueous and nonaqueous solutions based on accurate polarized partial charges. J. Chem. Theory Comput. 2007, 3, 2011−2033.

3369


Computing pKa Values in Different Solvents by ... - ACS Publications

Recommend Documents