We have developed a two-layer Molecules-in-Molecules (MIM2) fragmentation-based quantum chemical method including an efficient solvation model for the prediction of NMR chemical shifts with a target accuracy of ∼0.30 ppm for 1H and ∼2–3 ppm for 13C.
The Molecules-in-Molecules (MIM) fragmentation-based approach has been successfully used in previous studies to obtain the energies, optimized geometries, and spectroscopic properties of large molecular systems. The present work delineates a protocol to study the potential energy profiles for multistep chemical reactions using the MIM methodology. In a complex multistep chemical reaction, the fragmentation scheme needs to be changed as the reacting species transition into a new reaction step, resulting in a discontinuity in the potential energy curve of the reaction. In our approach, the fragmentation scheme for a particular step in a reaction is chosen on the basis of the nature of the bonding changes associated with that step. Thus, the reactant, transition state, and product are treated consistently throughout the reaction step, leading to an accurate energy barrier for that step. The discontinuity now occurs in describing the energies of reaction intermediates at the transition point between two reaction steps that are treated by two different fragmentation schemes. To address this issue, we propose a systematic procedure for obtaining continuous potential energy curves that are least shifted from their initial positions. The corrected MIM potential energy curves are continuous with activation energies preserved. Following this approach, energy profiles of complex reactions involving large molecular species can be obtained at high levels of theory with a reasonable computational cost.
Nitrosamines are in the cohort of concern (CoC) as determined by regulatory guidance. CoC compounds are considered highly potent carcinogens that need to be limited below the threshold of toxicological concern, 1.5 μg/day. Nitrosamines like NDMA and NDEA require strict control, while novel nitrosamine drug substance-related impurities (NDSRIs) may or may not be characterized as potent carcinogens. A risk assessment based on the structural features of NDSRIs is important in order to predict potency because they lack substance-specific carcinogenicity. Herein, we present a quantum mechanical (QM)-based analysis on structurally diverse sets of nitrosamines to better understand how structure influences the reactivity that could result in carcinogenicity. We describe the potency trend through activation energies corresponding to α-hydroxylation, aldehyde formation, diazonium intermediate formation, reaction with DNA base, and hydrolysis reactions, and other probable metabolic pathways associated with the carcinogenicity of nitrosamines. We evaluated activation energies for selected cases such as N-nitroso pyrrolidines, N-nitroso piperidines, N-nitroso piperazines, N-nitroso morpholines, N-nitroso thiomorpholine, N-methyl nitroso aromatic, fluorine-substituted nitrosamines, and substituted aliphatic nitrosamines. We compare these results to the recent framework of the carcinogenic potency characterization approach (CPCA) proposed by health authorities which is meant to give guidance on acceptable intakes (AI) for NDSRIs lacking substance-specific carcinogenicity data. We show examples where QM modeling and CPCA are aligned and examples where CPCA both underestimates and overestimates the AI. In cases where CPCA predicts high potency for NDSRIs, QM modeling can help better estimate an AI. Our results suggest that a combined mechanistic understanding of α-hydroxylation, aldehyde formation, hydrolysis, and reaction with DNA bases could help identify the structural features that underpin the potency of nitrosamines. We anticipate this work will be a valuable addition to the CPCA and provide a more analytical way to estimate AI for novel NDSRIs.
We have recently significantly expanded the applicability of our Molecules-in-Molecules (MIM) fragmentation method to large proteins by developing a three-layer model (MIM3) in which an accurate quantum-mechanical method is used in conjunction with a cost-effective, dispersion-corrected semiempirical model to overcome previous computational bottlenecks. In this work, we develop MIM3 as a structure-based drug design tool by application of the methodology for the accurate calculation of protein–ligand interaction energies. A systematic protocol is derived for the determination of the geometries of the protein–ligand complexes and to calculate their accurate interaction energies in the gas phase using MIM3. We also derive a simple and affordable procedure based on implicit solvation models and the ligand solvent-accessible surface area to approximate the ligand desolvation penalty in gas-phase interaction energy calculations. We have carefully assessed how closely such interaction energies, which are based on a single protein–ligand conformation, display correlations with the experimentally determined binding affinities. The performance of MIM3 was evaluated on a total of seven data sets comprising 89 protein–ligand complexes, all with experimentally known binding affinities, using a binding pocket involving a quantum region ranging in size from 250 to 600 atoms. The dispersion-corrected B97-D3BJ density functional, previously known to perform accurately for calculations involving non-covalent interactions, was used as the target level of theory for this work, with dispersion-corrected PM6-D3 as the semiempirical low level to incorporate the long-range interactions. Comparing directly to the experimental binding potencies, we obtain impressive correlations over all seven test sets, with an R2 range of 0.74–0.93 and a Spearman rank correlation coefficient (ρ) range of 0.83–0.93. Our results suggest that protein–ligand interaction energies are useful in predicting binding potency trends and validate the potential of MIM3 as a quantum-chemical structure-based drug design tool.
Experimental studies by Yamanouchi and co-workers indicate that an intense 40 fs 800 nm laser pulse can cause CH3OH+ to isomerizes during the pulse. The potential energy surfaces of methanol neutral, monocation, and singlet and triplet dication were explored using the CBS-APNO, CBS-QB3, CAM-B3LYP, and B3LYP levels of theory. Ab initio classical trajectories were calculated in the presence of a 2.9 × 1014 W/cm2 800 nm laser field for methanol monocation on the ground state potential energy surface using the CAM-B3LYP/6-31G(d,p) level of theory. With only zero point energy, CH3OH+ gained less than 15 kcal/mol from the 40 fs laser pulse, which was not enough to overcome any of the barriers for isomerization or fragmentation. To simulate extra energy deposited during the ionization process, 75, 100, and 125 kcal/mol of vibrational energy was added to the initial structures. After 400 fs, the distribution of product was CH2OH+ + H (79–81%), HCOH+ + H2 (9–13%), CH2OH2+ (1–3%), CH3+ + OH (1–3%), and CH2+ + H2O (<0.5%). The estimated kinetic energy releases are in accord with experimental findings. Experimental results using a probe pulse to ionize CH3OH+ to the dication showed substantial fraction C–O dissociation in both CH3OH+ and CH2OH2+ after the pulse. Because very few CH2OH2+ → CH2+ + H2O trajectories were seen in the simulation, the calculations suggest that some of the processes observed experimentally must occur on excited state surfaces or may be due to coupled nuclear-electron dynamics during the pump pulse.
Oxidative damage to DNA can lead to DNA–protein cross-links which can interfere with DNA transcription, replication, and repair. In experimental studies modeling oxidative damage to DNA, oxidation of guanosine by sulfate radical anion in the presence of lysine produced a mixture of lysine (Lys)-substituted spiroiminodihydantoins (Sp): ∼65% 5-Lys-Sp, ∼30% 8-Lys-Sp, and ∼5% 5,8-diLys-Sp. Pathways for formation of the lysine adducts during the oxidation of guanine by sulfate radical anions have been mapped out using B3LYP density functional theory and the SMD solvation model. Methylamine was used as a model for lysine, and imidazole served as a proton acceptor. The lowest barrier for methylamine reaction with guanine radical is addition at C8, yielding mainly 8-NHR-Sp and some 5,8-diNR-Sp. This is in good agreement with the cross-link ratios for mild oxidations mediated by type I photosensitizers such as benzophenone, but this is not in agreement with the product ratios for strong oxidants such as sulfate radical anion. The calculations explored pathways for oxidation of guanine by sulfate radical anion that produced guanine radical and radical cation and doubly oxidized guanine (Gox) and its cation. Sulfate radical anion can also oxidize methylamine to produce neutral methylamine radical (CH3NH•) after deprotonation. The calculations qualitatively reproduced the observed product ratio at pH 7 via a pathway involving the barrierless addition of methylamine radical at C5 and C8 of guanine radical. After C5 addition of methylamine radical, the lowest barrier is for H2O addition at C8 leading exclusively to 5-NHR-Sp. After C8 addition of methylamine radical, H2O and methylamine addition to C5 lead to 8-NHR-Sp and some 5,8-diNR-Sp.
A practical method for calculating the pKa's of selenols in aqueous solution has been developed by using density functional theory with the SMD solvation model and up to three explicit water molecules. The pKa's of 30 different organoselenols, 16 with known experimental pKa's, have been calculated by using three different functionals (ωB97XD, B3LYP, and M06-2X) and two basis sets (6-31+G(d,p) and 6-311++G(d,p)). Calculations using ωB97XD and B3LYP with SMD solvation without explicit waters are found to have errors of 3-6 pKa units; the errors with M06-2X are somewhat smaller. One explicit water interacting with the selenium reduces the calculated pKa's by 1-2 pKa units along with improving the slope and intercept of the fit of the calculated pKa's to experiment. The best results for SMD/M06-2X/6-31+G(d,p) are with one explicit water (MSE = -0.08 ± 0.37 and MUE = 0.32 ± 0.37) whereas ωB97XD and B3LYP still have errors larger than 2 pKa units. The best results for ωB97XD and B3LYP with 6-31+G(d,p) are obtained by using three explicit waters (MSE = 0.36 ± 0.24 and 0.34 ± 0.25, respectively) and a fit to experiment yields a slope of 1.06 with a zero intercept. The errors for M06-2X/6-31+G(d,p) with three explicit waters are significantly larger (-3.59 ± 0.45) because it overstabilizes the anions. On the basis of the molecules studied here, M06-2X/6-31+G(d,p) with one explicit water and ωB97XD/6-31+G(d,p) and B3LYP/6-31+G(d,p) with three explicit waters and the SMD solvation model can produce reliable pKa's for the substituted selenols.
 scale(5.56106)'%3e%3ccircle r='1'/%3e%3cg id='d'%3e%3cg id='c'%3e%3cpath id='b' d='M.195-.98 0-1.39l-.195.408' transform='rotate(11.25)'/%3e%3cuse xlink:href='%23b' transform='rotate(22.5)'/%3e%3cuse xlink:href='%23b' transform='rotate(45)'/%3e%3c/g%3e%3cuse xlink:href='%23c' transform='rotate(67.5)'/%3e%3c/g%3e%3cuse xlink:href='%23d' transform='scale(-1 1)'/%3e%3c/g%3e%3cg transform='translate(0 58.787) scale(8.1434)'%3e%3ccircle r='1'/%3e%3cg id='g'%3e%3cg id='f'%3e%3cpath id='e' d='M.259.966 0 1.576l-.259-.61'/%3e%3cuse xlink:href='%23e' transform='rotate(180)'/%3e%3c/g%3e%3cuse xlink:href='%23f' transform='rotate(90)'/%3e%3c/g%3e%3cuse xlink:href='%23g' transform='rotate(30)'/%3e%3cuse xlink:href='%23g' transform='rotate(60)'/%3e%3c/g%3e%3c/g%3e%3c/svg%3e)
'%3e%3ccircle r='20' fill='%23008'/%3e%3ccircle r='17.5' fill='%23fff'/%3e%3ccircle r='3.5' fill='%23008'/%3e%3cg id='d'%3e%3cg id='c'%3e%3cg id='b'%3e%3cg id='a' fill='%23008'%3e%3ccircle r='.875' transform='rotate(7.5 -8.75 133.5)'/%3e%3cpath d='M0 17.5.6 7 0 2l-.6 5L0 17.5z'/%3e%3c/g%3e%3cuse xlink:href='%23a' transform='rotate(15)'/%3e%3c/g%3e%3cuse xlink:href='%23b' transform='rotate(30)'/%3e%3c/g%3e%3cuse xlink:href='%23c' transform='rotate(60)'/%3e%3c/g%3e%3cuse xlink:href='%23d' transform='rotate(120)'/%3e%3cuse xlink:href='%23d' transform='rotate(-120)'/%3e%3c/g%3e%3c/svg%3e)