Prediction of acidity in acetonitrile solution with COSMO-RS

The COSMO-RS method, a combination of the quantum chemical dielectric continuum solvation model COSMO with a statistical thermodynamics treatment for realistic solvation simulations, has been used for the prediction of pK(a) values in acetonitrile. For a variety of 93 organic acids, the directly calculated values of the free energies of dissociation in acetonitrile showed a very good correlation with the pK(a) values (r(2) = 0.97) in acetonitrile, corresponding to a standard deviation of 1.38 pK(a) units. Thus, we have a prediction method for acetonitrile pK(a) with the intercept and the slope as the only adjusted parameters. Furthermore, the pK(a) values of CH acids yielding large anions with delocalized charge can be predicted with a rmse of 1.12 pK(a) units using the theoretical values of slope and intercept resulting in truly ab initio pK(a) prediction. In contrast to our previous findings on aqueous acidity predictions the slope of the experimental pK(a) versus theoretical DeltaG(diss) was found to match the theoretical value 1/RT ln (10) very well. The predictivity of the presented method is general and is not restricted to certain compound classes. However, a systematic correction of -7.5 kcal mol(-1) is required for compounds that do not allow electron-delocalization in the dissociated anion. The prediction model was tested on a diverse test set of 129 complex multifunctional compounds from various sources, reaching a root mean square deviation of 2.10 pK(a) units.     

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?

A popular method for predicting pK(a) values for organic molecules in aqueous solution is to establish empirical linear least-squares fits between calculated deprotonation energies and known experimental pK(a) values. In virtually all such calculations, the empirically observed slope of the pK(a) vs. ΔE fit is significantly less than the theoretical value, 1/(2.303RT) (which is 0.73 mol/kcal at room temperature). In our own continuum solvation calculations (Zhang et al., J Phys Chem A 2010, 114, 432), the empirical slope for carboxylic acids was only 0.23 mol/kcal, despite the excellent fit to the experimental pK(a) values. There has been much speculation about the reason for this phenomenon. Although the ΔE - pK(a) relation neglects entropic effects, these are expected to largely cancel. The most likely cause for the strange behavior of the fitted slope is explicit solute-solvent (water) interactions, especially involving the ions, which cannot be described accurately by continuum solvation models. We used our previously developed pK(a) protocol (OLYP/6-311+G(d,p)//3-21G(d) with the COSMO solvation model) to investigate the effect of adding one or two explicit water molecules to the system. The slopes for organic acids (especially carboxylic acids) are much closer to the theoretical value when explicit water molecules are added to both the neutral molecule and the anion. However, explicit water molecules have almost no effect on the slopes for organic bases. Adding explicit water molecules to the ions only produces intermediate results. Unfortunately, linear fits involving explicit water molecules have much larger errors than with continuum solvation models alone and are also much more expensive. Consequently, they are not suitable for large-scale pK(a) calculations. The results compared with literature values showed that our predicted pK(a) s are more accurate.     

Calculations of pKa of superacids in 1,2-dichloroethane

Acidity calculations for some CH and NH superacids in 1,2-dichloroethane (DCE) were carried out using SMD and COSMO-RS continuum solvation models. After comparing the results of calculations with respective experimental pK(a) values it was found that the performance of SMD/M05-2X/6-31G* method is characterized by the mean unsigned error (MUE) of 0.5 pK(a) units and the slope of regression line of 0.915. The similar SMD/B3LYP/6-31G* approach was slightly less successful. The strong correlation over entire data set is confirmed by R(2) values of 0.990 and 0.984 for M05-2X and B3LYP functionals, respectively. The COSMO-RS method, while providing the value of the linear regression line slope similar to the corresponding values from SMD approach, characterized by rather loose correlation (R(2) = 0.823, MUE = 1.7 pK(a) units) between calculated and experimental pK(a) values in DCE solution.     

Dependence of pKa on solute cavity for diprotic and triprotic acids

A systematic study of ΔG(aq)/pK(a) for monoprotic, diprotic, and triprotic acids has been carried out based on DFT/aug-cc-pVTZ combined with CPCM and SMD solvation modeling. All DFT/cavity set combinations considered showed similar accuracy for ΔG(aq)(1)/pK(a1) (70% within ±2.5 kcal mol(-1) of experiment) while only the M05-2X/Pauling cavity combination gave reasonable results for ΔG(aq)(2)/pK(a2) when both pK(a) values are separated by more than three units (70% within ±5.0 kcal mol(-1) of experiment). The choice of experimental data is critical to the interpretation of the calculated accuracy especially for several inorganic acids. For the calculation of ΔG(aq)(3)/pK(a3), the larger experimental uncertainty and an unrealistic orbital population of diffuse function for trianions in the gas phase hinders an evaluation of the predictive performance. We find the M05-2X functional with the Pauling cavity set is the best choice for ΔG(aq)(2)/pK(a2) prediction in aqueous media while all DFT/cavity sets considered were competitive for ΔG(aq)(1)/pK(a1).     

Prediction of the dissociation constant pKa of organic acids from local molecular parameters of their electronic ground state

A quantum chemical method has been developed to estimate the dissociation constant pK(a) of organic acids from their neutral molecular structures by employing electronic structure properties. The data set covers 219 phenols (including 29 phenols with intramolecular H-bonding), 150 aromatic carboxylic acids, 190 aliphatic carboxylic acids, and 138 alcohols, with pK(a) varying by 16 units (0.38-16.80). Optimized ground-state geometries employing the semiempirical AM1 Hamiltonian have been used to quantify the site-specific molecular readiness to donate or accept electron charge in terms of both charge-associated energies and energy-associated charges, augmented by an ortho substitution indicator for aromatic compounds. The resultant regression models yield squared correlation coefficients (r(2)) from 0.82 to 0.90 and root-mean-square errors (rms) from 0.39 to 0.70 pK(a) units, corresponding to an overall (subset-weighted) r(2) of 0.86. Simulated external validation, leave-10%-out cross-validation and target value scrambling demonstrate the statistical robustness and prediction power of the derived model suite. The low intercorrelation with prediction errors from the commercial ACD package provides opportunity for a consensus model approach, offering a pragmatic way for further increasing the confidence in prediction significantly. Interestingly, inclusion of calculated free energies of aqueous solvation does not improve the prediction performance, probably because of the limited precision provided by available continuum-solvation models.     

Ab Initio Quantum Mechanical Study on the Origin of the pK(a) Differences of the Proton Sponges 1,8-Bis(dimethylamino)naphthalene, 1,8-Bis(dimethylamino)-2,7-dimethoxynaphthalene, 1,6-Dimethyl-1,6-diazacyclodecane, and 1,6-Diazabicyclo[4.4.4]tetradecane

Ab initio quantum mechanical calculations were used in studying the origin of the exceptionally high basicities of four diamines (13-16) with pK(a1) values ranging from 12.1 to 25. The computational approach involved the calculation of the gas-phase proton affinities of the molecules studied at the MP2/6-31G//HF/6-31G level and the solvation energies with the polarizable continuum model at the HF/6-31G level. The calculated gas-phase and aqueous-phase proton affinities of a structurally diverse series of amines were compared with the corresponding experimental gas-phase proton affinities and pK(a1) values. The calculated values were found to be in reasonable agreement with the experimental ones. The high basicities of the studied diamines were found to originate from the nitrogen lone-pair repulsion, solvation effects, and strong intramolecular hydrogen bonds. Each of these factors were found to be able to increase the pK(a1) values of the high-basicity diamines by 2-6 pK(a) units. The relative contributions of the factors varied between the compounds. The nitrogen lone-pair repulsion was estimated to be the most important factor in increasing the pK(a1) values. In addition, barriers for proton transfers between the nitrogens of selected diamines were calculated, and comparison was made between the barrier heights and the geometries of the diamines.     

Protonation studies of modified adenine and adenine nucleotides by theoretical calculations and (15)N NMR

The acid/base character of nucleobases affects phenomena such as self-association, interaction with metal ions, molecular recognition by proteins, and nucleic acid base-pairing. Therefore, the investigation of proton-transfer equilibria of natural and synthetic nucleos(t)ides is of great importance to obtain a deeper understanding of these phenomena. For this purpose, a set of ATP prototypes was investigated using (15)N NMR spectroscopy, and the corresponding adenine bases were investigated by theoretical calculations. (15)N NMR measurements provided not only acidity constants but also information on the protonation site(s) on the adenine ring and regarding the ratio of the singly protonated species in equilibrium. Substituents of different nature and position on the adenine ring did not change the preferred protonation site, which remained N1. However, for 2-thioether-ATP derivatives a mixed population of N1 and N7 singly protonated species was observed. Reduction of basicity of 0.4-1 pK(a) units relative to ATP was also observed for all evaluated ATP derivatives, except for 2-Cl-ATP, for which K(a) was ca. 10,000-fold lower. To explain the substitution-dependent variations in the experimental pK(a) values of the ATP analogues, gas-phase proton affinities (PA), Delta Delta G(hyd), and pK(a) values of the corresponding adenine bases were calculated using quantum mechanical methods. The computed PA and Delta Delta G(hyd) values successfully explained the experimental pK(a) values. A computational procedure for the prediction of accurate pK(a) values was developed using density functional theory and polarizable continuum model calculations. In this procedure, we developed a set of parameters for the polarizable continuum model that was fitted to reproduce experimental pK(a) values of nitrogen heterocycles. This method is proposed for the prediction of pK(a) values and protonation site(s) of purine analogues that have not been synthesized or analyzed.     

Influence of water content on the acidities in acetonitrile. Quantifying charge delocalization in anions

The effect of traces of water on the relative strengths of acids (ΔpK(a) values) in acetonitrile was quantitatively evaluated experimentally and computationally (COSMO-RS). Water affects first of all the anions by selective solvation. Expectedly, the more localized is the charge in acid anions the higher is the effect of water. The energetic effect of increasing water content from 0 to ca. 10,000 ppm on solvation enthalpies of anions ranged from 0.2-0.4 kcal mol?1 (anions with delocalized charges) to 15 kcal mol?1 in the case of the highly charge-localized acetate ion. In the case of ΔpK(a) values the change ranges from 0.01 to ca. 1.7 pK(a) units (acid pair involving acetic acid). The COSMO-RS method was found to satisfactorily describe the trends in ΔpK(a) values. To quantify the extent of charge localization/delocalization in anions a parameter, weighted average positive σ (WAPS), was introduced, which can be conveniently computed using the COSMO approach. WAPS characterizes the distribution of charge density across the molecular surface and was found to correlate well with the extent of water influence on the dissociation of the respective acid.     

Accurate pK(a) calculations for carboxylic acids using complete basis set and Gaussian-n models combined with CPCM continuum solvation methods

Complete Basis Set and Gaussian-n methods were combined with CPCM continuum solvation methods to calculate pK(a) values for six carboxylic acids. An experimental value of -264.61 kcal/mol for the free energy of solvation of H(+), DeltaG(s)(H(+)), was combined with a value for G(gas)(H(+)) of -6.28 kcal/mol to calculate pK(a) values with Cycle 1. The Complete Basis Set gas-phase methods used to calculate gas-phase free energies are very accurate, with mean unsigned errors of 0.3 kcal/mol and standard deviations of 0.4 kcal/mol. The CPCM solvation calculations used to calculate condensed-phase free energies are slightly less accurate than the gas-phase models, and the best method has a mean unsigned error and standard deviation of 0.4 and 0.5 kcal/mol, respectively. The use of Cycle 1 and the Complete Basis Set models combined with the CPCM solvation methods yielded pK(a) values accurate to less than half a pK(a) unit.     

Theoretical evaluation of pK(a) in phosphoranes: implications for phosphate ester hydrolysis

Knowledge of the pK(a) of phosphoranes is important for the interpretation of phosphate ester hydrolysis. Calculated pK(a)'s of the model phosphorane, ethylene phosphorane, are reported. The method of calculation is based on the use of dimethyl phosphate as a reference state for evaluating relative pK(a) values, and on the optimization of the oxygen and acidic hydrogen van der Waals radii to give reasonable pK(1)(a), pK(2)(a), and pK(3)(a) for phosphoric acid in solution. Density functional theory is employed to calculate the gas-phase protonation energies, and continuum dielectric methods are used to determine the solvation corrections. The calculated pK(1)(a) and p(2)(a) for the model phosphorane are 7.9 and 14.3, respectively. These values are within the range of proposed experimental values, 6.5-11.0 for pK(1)(a), and 11.3-15.0 for pK(2)(a). The mechanistic implications of the calculated pK(a)'s are discussed.     

pKa prediction from an ab initio bond length: part 3--benzoic acids and anilines

The prediction of pK(a) from a single ab initio bond length has been extended to provide equations for benzoic acids and anilines. The HF/6-31G(d) level of theory is used for all geometry optimisations. Similarly to phenols (Part 2 of this series of publications), the meta-/para-substituted benzoic acids can be predicted from a single model constructed from one bond length. This model had an impressive RMSEP of 0.13 pK(a) units. The prediction of ortho-substituted benzoic acids required the identification of high-correlation subsets, where the compounds in the same subset have at least one of the same (e.g. halogens, hydroxy) ortho substituent. Two pK(a) equations are provided for o-halogen benzoic acids and o-hydroxybenzoic acids, where the RMSEP values are 0.19 and 0.15 pK(a) units, respectively. Interestingly, the bond length that provided the best model differed between these two high-correlation subsets. This demonstrates the importance of investigating the most predictive bond length, which is not necessarily the bond involving the acid hydrogen. Three high-correlation subsets were identified for the ortho-substituted anilines. These were o-halogen, o-nitro and o-alkyl-substituted aniline high-correlation subsets, where the RMSEP ranged from 0.23 to 0.44 pK(a) units. The RMSEP for the meta-/para-substituted aniline model was 0.54 pK(a) units. This value exceeded our threshold of 0.50 pK(a) units and was higher than both the m-/p-benzoic acids in this work and the m-/p-phenols (RMSEP = 0.43) of Part 2. Constructing two separate models for the meta- and para- substituted anilines, where RMSEP values of 0.63 and 0.33 pK(a) units were obtained respectively, revealed it was the meta-substituted anilines that caused the large RMSEP value. For unknown reasons the RMSEP value increased with the addition of a further twenty meta-substituted anilines to this model. The C-N bond always produced the best correlations with pK(a) for all the high-correlation subsets. A higher level of theory and an ammonia probe improved the statistics only marginally for the hydroxybenzoic acid high-correlation subsets.     

First-principles prediction of the pK(a)s of anti-inflammatory oxicams

The gas- and aqueous-phase acidities of a series of oxicams have been computed by combining M05-2X/6-311+G(3df,2p) gas-phase free energies with solvation free energies from the CPCM-UAKS, COSMO-RS, and SMD solvent models. To facilitate accurate gas-phase calculations, a benchmarking study was further carried out to assess the performance of various density functional theory methods against the high-level composite method G3MP2(+). Oxicams are typically diprotic acids, and several tautomers are possible in each protonation state. The direct thermodynamic cycle and the proton exchange scheme have been employed to compute the microscopic pK(a)s on both solution- and gas-phase equilibrium conformers, and these were combined to yield the macroscopic pK(a) values. Using the direct cycle of pK(a) calculation, the CPCM-UAKS model delivered reasonably accurate results with MAD ~ 1, whereas the SMD and COSMO-RS models' performance was less satisfactory with MAD ~ 3. Comparison with experiment also indicates that direct cycle calculations based on solution conformers generally deliver better accuracy. The proton exchange cycle affords further improvement for all solvent models through systematic error cancellation and therefore provides better reliability for the pK(a) prediction of compounds of these types. The latter approach has been applied to predict the pK(a)s of several recently synthesized oxicam derivatives.     

Absolute pK(a) determinations for substituted phenols

The CBS-QB3 method was used to calculate the gas-phase free energy difference between 20 phenols and their respective anions, and the CPCM continuum solvation method was applied to calculate the free energy differences of solvation for the phenols and their anions. The CPCM solvation calculations were performed on both gas-phase and solvent-phase optimized structures. Absolute pK(a) calculations with solvated phase optimized structures for the CPCM calculations yielded standard deviations and root-mean-square errors of less than 0.4 pK(a) unit. This study is the most accurate absolute determination of the pK(a) values of phenols, and is among the most accurate of any such calculations for any group of compounds. The ability to make accurate predictions of pK(a) values using a coherent, well-defined approach, without external approximations or fitting to experimental data, is of general importance to the chemical community. The solvated phase optimized structures of the anions are absolutely critical to obtain this level of accuracy, and yield a more realistic charge separation between the negatively charged oxygen and the ring system of the phenoxide anions.     

On the role of solvent in acid-base equilibria of diuretics in acetonitrile-water mixed solvents

The dissociation constants of representative loop (furosemide), thiazide (chlorthiazide and trichlormethiazide) and potassium sparing (amiloride) diuretics in 10, 30, 40, 50 and 70% (w/w) acetonitrile-water mixtures at 298.15 K were determined, in accordance with IUPAC procedures. The variation in pK(a) values over the whole composition range can be explained by preferential solvation and the structural features in acetonitrile-water mixtures. Correlations between pK(a) values and various bulk and solvatochromic properties of the solvents were calculated. The linear solvation energy relationship (LSER) method was applied. The resulting equations allowed us to calculate pK(a) values for the diuretics studied in any acetonitrile-water mixture up to 70% (w/w) and should help to clarify the acid-base behaviour of diuretics in the widely-used acetonitrile-water mixed solvents.     

Protonation equilibria of quinolone antibacterials in acetonitrile-water mobile phases used in LC

Ionization constants of nine quinolone antibacterials in acetonitrile-water mixtures containing 0, 5.5, 10, 16.3, 25, 30, 40, 50 and 70% (w/w) acetonitrile were obtained and assignment of these pK values to the several potentially ionizable functional groups was made. The variation of the pK values obtained over the whole composition range studied can be explained by consideration of the preferential solvation of electrolytes in acetonitrile-water mixtures. In order to obtain pK values in any of the unlimited number of possible binary solvent acetonitrile-water mixtures, relationships between pK values and different bulk properties (such as dielectric constant) were examined. The linear solvation energy relationships method, LSER, was applied to study the correlation of pK values with the solvatochromic parameters of acetonitrile-water mixtures. The equations obtained allow calculation of the pK values of the quinolone antimicrobials in any acetonitrile-water mixtures up to 70% (w/w) and thus permit the knowledge of the acid-base behaviour of these important antimicrobials in the widely used acetonitrile-water media.     

Theoretical prediction of pKa values of seleninic, selenenic, sulfinic, and carboxylic acids by quantum-chemical methods

Aqueous acid dissociation constants of substituted areneseleninic, areneselenenic, arenesulfinic, and benzoic acids are calculated by ab initio (MP2) and DFT (B3LYP) methods in combination with bulk solvation models (IEFPCM, CRSrad) from appropriate thermodynamic cycles. Mean absolute deviations (MAD) between experimental and calculated pK(a) values are quite large for basis sets without diffuse functions; however, trends are reasonably well described. Best agreement with experiment as described by MAD as well as correlation coefficient and slope of the correlation equation pK(a) = a*ΔG(calc)/RT ln(10) + b is obtained with the CPCM solvation model using the defaults optimized within COSMO-RS (CRSrad; MAD = 1.54, R(2) = 0.94, a = 0.83). Sulfenic (selenenic) acid tautomers are significantly more stable than the corresponding sulfoxide (selenoxide) forms.     

Acid-base equilibria in gamma-butyrolactone studied by use of pH-ISFETs

pH-ISFETs were used in the study of acid-base equilibria in gamma-butyrolactone (GBL). After the spectrophotometric determination of the pK(a) value of 3,5-dichloropicric acid, the pK(a) values and homo-conjugation constants of various acids (including the conjugate acids of bases) were determined potentiometrically using a Ta(2)O(5)-type pH-ISFET. The values of pK(a) in GBL were in a linear relation with those in propylene carbonate (PC) and 1.0 units smaller on average. The difference in pK(a) between GBL and PC was mainly attributable to the difference in proton solvation. The autoprotolysis constant of GBL, roughly estimated by a rapid titration with a Si(3)N(4)-ISFET, was about 30 on the pK(SH) scale. A comparative study was made of the response speeds of the Ta(2)O(5)- and Si(3)N(4)-type pH-ISFETs and a conventional pH-glass electrode. The result was Si(3)N(4)-ISFET > Ta(2)O(5)-ISFET > glass electrode. Because GBL is not stable against acids and bases, the use of pH-ISFETs was much more convenient than the use of the conventional glass electrode.     

Determination of the overlapping pK(a) values of chrysin using UV-vis spectroscopy and ab initio methods

The overlapping pK(a) values of 5,7-dihydroxyflavone (chrysin) in EtOH-water solutions were determined by means of a UV-vis spectroscopic method that uses absorbance diagrams, at constant ionic strength (0.050 M) and temperature (25.0+/-0.1 degrees C). It was observed that the pK(a) values increase when the polarity-polarizability and solvation abilities of the reaction medium decrease. In order to calculate the pK(a1) and pK(a2) of chrysin in pure water, various relationships between the determined pK(a) and properties of solvents (relative permittivity, alpha-parameter of Taft and parameter Acity), are proposed. Moreover, with the aim of explaining the first pK(a1) value obtained, the molecular conformations and solute-solvent interactions of the 7(O(-))chrysinate monoanion were also investigated, using ab initio methods. Several ionization reactions and equilibria in water, which possesses a high hydrogen-bond-donor ability, are proposed. These reactions and equilibria constituted the necessary theoretical basis to calculate the first acidity constant of chrysin. The HF/6-31G(d) and HF/6-31+G(d) methods were used for calculations. Tomasi's method was used to analyze the formation of intermolecular hydrogen bonds between the 7(O(-))chrysinate monoanion and water molecules. It was proposed that in alkaline aqueous solutions the monoanion of chrysin is solvated with one water molecule. The agreement between the experimental and theoretical pK(a1) values provides good support for the acid-base reactions proposed in this paper.     

pKa of acetate in water: a computational study

Several computational methods including the conductor-like polarizable continuum model, CPCM with both UAKS and UAHF cavities, Cramer and Truhlar's generalized Born solvation model, SM5.4(AM1), SM5.4(PM3), and SM5.43R(mPW1PW91/6-31+G(d)), and mixed QM/MM-Ewald simulations were used to calculate the pK(a) values of acetate and bicarbonate anions in aqueous solution. This work provided a critical and comprehensive assessment of the quality of these theoretical models in the calculation of aqueous solvation free energies for the singly charged acetate and bicarbonate ions, as well as the doubly charged acetate dianion and carbonate dianion. It was shown that QM/MM-Ewald simulations could give an accurate and consistent evaluation of the pK(a) values of acetate and bicarbonate based on both the relative and absolute pK(a) formulas, while other methods could yield satisfactory results only for certain calculations. However, this does not mean that the current QM/MM-Ewald protocol is superior to other methods. The useful information obtained in this investigation is that both the absolute and relative pK(a) formulas should better be tested in accurate calculations of pK(a) values based on any methods.     

FAMBE-pH: a fast and accurate method to compute the total solvation free energies of proteins

A fast and accurate method to compute the total solvation free energies of proteins as a function of pH is presented. The method makes use of a combination of approaches, some of which have already appeared in the literature; (i) the Poisson equation is solved with an optimized fast adaptive multigrid boundary element (FAMBE) method; (ii) the electrostatic free energies of the ionizable sites are calculated for their neutral and charged states by using a detailed model of atomic charges; (iii) a set of optimal atomic radii is used to define a precise dielectric surface interface; (iv) a multilevel adaptive tessellation of this dielectric surface interface is achieved by using multisized boundary elements; and (v) 1:1 salt effects are included. The equilibrium proton binding/release is calculated with the Tanford-Schellman integral if the proteins contain more than approximately 20-25 ionizable groups; for a smaller number of ionizable groups, the ionization partition function is calculated directly. The FAMBE method is tested as a function of pH (FAMBE-pH) with three proteins, namely, bovine pancreatic trypsin inhibitor (BPTI), hen egg white lysozyme (HEWL), and bovine pancreatic ribonuclease A (RNaseA). The results are (a) the FAMBE-pH method reproduces the observed pK a's of the ionizable groups of these proteins within an average absolute value of 0.4 p K units and a maximum error of 1.2 p K units and (b) comparison of the calculated total pH-dependent solvation free energy for BPTI, between the exact calculation of the ionization partition function and the Tanford-Schellman integral method, shows agreement within 1.2 kcal/mol. These results indicate that calculation of total solvation free energies with the FAMBE-pH method can provide an accurate prediction of protein conformational stability at a given fixed pH and, if coupled with molecular mechanics or molecular dynamics methods, can also be used for more realistic studies of protein folding, unfolding, and dynamics, as a function of pH.