Stefan Balaz

Multimode ligand binding in receptor site modeling: Implementation in CoMFA

Receptor site modeling methods usually use one binding mode (conformation and/or orientation) for each ligand in a 1:1 complex with receptor. Multiple modes should be considered instead because (1). they have frequently been observed experimentally; (2). in a series, ligands can bind in single yet different modes; and (3). a series may only exhibit one but unknown mode and a few plausible modes must be examined. For multimode binding, the observed ligand/receptor association constant is the sum of the association constants that characterize individual binding modes. This relation, when applied to Comparative Molecular Field Analysis (CoMFA), results in a dependence of the observed binding energy on the probe energies that is nonlinear in optimized parameters. The dependence was linearized to allow parameter optimization by the partial least-squares method that was used iteratively until self-consistency. In addition to the standard CoMFA output, the procedure objectively selects one or a few optimal binding modes out of a dozen or more modes that are considered for each ligand. The approach was applied to published data for binding of 34 polychlorinated dibenzofurans to the aryl hydrocarbon receptor. Descriptive and predictive abilities of the 16-mode model were significantly better than for the one-, two-, and four-mode models. Predominantly, edge-aligned modes were selected that are seldom used in CoMFA. Since inclusion of multimode binding only changes the form of the correlation equation and does not affect the number of optimized parameters, the improvement is believed to be due to a more realistic description.

QM/MM Linear Response method distinguishes ligand affinities for closely related metalloproteins

Design of selective ligands for closely related targets is becoming one of the most important tasks in the drug development. New tools, more precise than fast scoring functions and less demanding than sophisticated Free Energy Perturbation methods, are necessary to help accomplish this goal. The methods of intermediate complexity, characterizing individual contributions to the binding energy, have been an area of intense research in the past few years. Our recently developed quantum mechanical/molecular mechanical (QM/MM) modification of the Linear Response (LR) method describes the binding free energies as the sum of empirically weighted contributions of the QM/MM interaction energies and solvent‐accessible surface areas for the time‐averaged structures of hydrated complexes, obtained by molecular dynamics (MD) simulations. The method was applied to published data on 27 inhibitors of matrix metalloproteinase‐3 (MMP‐3). The two descriptors explained 90% of variance in the inhibition constants with RMSE of 0.245 log units. The QM/MM treatment is indispensable for characterization of the systems lacking suitable force‐field expressions. In this case, it provided characteristics of H‐bonds of the inhibitors to Glu202, charges of binding site atoms, and accurate coordination geometries of the ligands to catalytic zinc. The geometries were constrained during the MD simulations, which characterized conformational flexibility of the complexes and helped in the elucidation of the binding differences for related compounds. A comparison of the presented QM/MM LR results with those previously published for inhibition of MMP‐9 by the same set of ligands showed that the QM/MM LR approach was able to distinguish subtle differences in binding affinities for MMP‐3 and MMP‐9, which did not exceed one order of magnitude. This precision level makes the approach a useful tool for design of selective ligands to similar targets, because the results can be safely extrapolated to maximize selectivity. Proteins 2007.

Binding affinity prediction for ligands and receptors forming tautomers and ionization species: inhibition of mitogen-activated protein kinase-activated protein kinase 2 (MK2).

Treatment of ionization and tautomerism of ligands and receptors is one of the unresolved issues in structure-based prediction of binding affinities. Our solution utilizes the thermodynamic master equation, expressing the experimentally observed association constant as the sum of products, each valid for a specific ligand–receptor species pair, consisting of the association microconstant and the fractions of the involved ligand and receptor species. The microconstants are characterized by structure-based simulations, which are run for individual species pairs. Here we incorporated the multispecies approach into the QM/MM linear response method and used it for structural correlation of published inhibition data on mitogen-activated protein kinase (MAPK)-activated protein kinase (MK2) by 66 benzothiophene and pyrrolopyridine analogues, forming up to five tautomers and seven ionization species under experimental conditions. Extensive cross-validation showed that the resulting models were stable and predictive. Inclusion of all tautomers and ionization ligand species was essential: the explained variance increased to 90% from 66% for the single-species model.

Subcellular pharmacokinetics and its potential for library focusing.

Subcellular pharmacokinetics (SP) optimizes biology-related factors in the design of libraries for high throughput screening by defining comparatively narrow ranges of properties (lipophilicity, amphiphilicity, acidity, reactivity, 3D-structural features) of the included compounds. The focusing ensures appropriate absorption, distribution, metabolism, excretion, and toxicity (ADMET) in those test biosystems, which are more complex than isolated receptors, and in humans. The SP deploys conceptual models that include transport and accumulation in a series of membranes, protein binding, hydrolysis, and other reactions with cell constituents. The kinetics of drug disposition is described as a non-linear disposition function of drug structure and properties. The SP capabilities are illustrated here using a model-based quantitative structure-activity relationship of toxicity of phenolic compounds against Tetrahymena pyriformis as dependent on lipophilicity and acidity. The resulting SP models clearly outperform empirical models in predictive ability outside the parameter space, as revealed by the leave-extremes-out cross-validation technique with omission of compounds beyond pre-defined lipophilicity and acidity ranges. The SP models do not change substantially if the parameters space is shrunk within some limits. In contrast, the shapes of empirical models vary widely depending upon the fraction of the data set used for their optimization. Once calibrated for a given biosystem, the SP models provide a detailed recipe for tailoring the drug properties to ensure optimum ADMET. The focusing is more accurate than with traditional empirical QSAR studies, assessment of drug-likeness, or the rules for identification of compounds with permeability problems.

Improved estimation of ligand–macromolecule binding affinities by linear response approach using a combination of multi-mode MD simulation and QM/MM methods

Structure-based predictions of binding affinities of ligands binding to proteins by coordination bonds with transition metals, covalent bonds, and bonds involving charge re-distributions are hindered by the absence of proper force fields. This shortcoming affects all methods which use force–field-based molecular simulation data on complex formation for affinity predictions. One of the most frequently used methods in this category is the Linear Response (LR) approach of Åquist, correlating binding affinities with van der Waals and electrostatic energies, as extended by Jorgensen’s inclusion of solvent-accessible surface areas. All these terms represent the differences, upon binding, in the ensemble averages of pertinent quantities, obtained from molecular dynamics (MD) or Monte Carlo simulations of the complex and of single components. Here we report a modification of the LR approach by: (1) the replacement of the two energy terms through the single-point QM/MM energy of the time-averaged complex structure from an MD simulation; and (2) a rigorous consideration of multiple modes (mm) of binding. The first extension alleviates the force-field related problems, while the second extension deals with the ligands exhibiting large-scale motions in the course of an MD simulation. The second modification results in the correlation equation that is nonlinear in optimized coefficients, but does not lead to an increase in the number of optimized coefficients. The application of the resulting mm QM/MM LR approach to the inhibition of zinc-dependent gelatinase B (matrix metalloproteinase 9) by 28 hydroxamate ligands indicates a significant improvement of descriptive and predictive abilities.

Drug-membrane interactions studied in phospholipid monolayers adsorbed on nonporous alkylated microspheres.

Characterization of interactions with phospholipids is an integral part of the in vitro profiling of drug candidates because of the roles the interactions play in tissue accumulation and passive diffusion. Currently used test systems may inadequately emulate the bilayer core solvation properties (immobilized artificial membranes [IAM]), suffer from potentially slow transport of some chemicals (liposomes in free or immobilized forms), and require a tedious separation (if used for free liposomes). Here the authors introduce a well-defined system overcoming these drawbacks: nonporous octadecylsilica particles coated with a self-assembled phospholipid monolayer. The coating mimics the structure of the headgroup region, as well as the thickness and properties of the hydrocarbon core, more closely than IAM. The monolayer has a similar transition temperature pattern as the corresponding bilayer. The particles can be separated by filtration or a mild centrifugation. The partitioning equilibria of 81 tested chemicals were dissected into the headgroup and core contributions, the latter using the alkane/water partition coefficients. The deconvolution allowed a successful prediction of the bilayer/water partition coefficients with the standard deviation of 0.26 log units. The plate-friendly assay is suitable for high-throughput profiling of drug candidates without sacrificing the quality of analysis or details of the drug-phospholipid interactions.

Structural determinants of drug partitioning in surrogates of phosphatidylcholine bilayer strata.

The knowledge of drug concentrations in bilayer headgroups, core, and at the interface between them is a prerequisite for quantitative modeling of drug interactions with many membrane-bound transporters, metabolizing enzymes and receptors, which have the binding sites located in the bilayer. This knowledge also helps understand the rates of trans-bilayer transport because balanced interactions of drugs with the bilayer strata lead to high rates, while excessive affinities for any stratum cause a slowdown. Experimental determination of bilayer location is so tedious and costly that the data are only available for some fifty compounds. To extrapolate these valuable results to more compounds at a higher throughput, surrogate phases have been used to obtain correlates of the drug affinities for individual strata. We introduced a novel system, consisting of a diacetyl phosphatidylcholine (DAcPC) solution with the water content of the fluid bilayer as the headgroup surrogate and n-hexadecane (C16) representing the core. The C16/DAcPC partition coefficients were measured for 113 selected compounds, containing structural fragments that are frequently occurring in approved drugs. The data were deconvoluted into the ClogP-based fragment solvation characteristics and processed using a solvatochromic correlation. Increased H-bond donor ability and excess molar refractivity of compounds promote solvation in the DAcPC phase as compared to bulk water, contrary to H-bond acceptor ability, dipolarity/polarizability, and volume. The results show that aromates have more balanced distribution in bilayer strata, and thus faster trans-bilayer transport, than similar alkanes. This observation is in accordance with the frequent occurrence of aromatic rings in approved drugs and with the role of rigidity of drug molecules in promoting intestinal absorption. Bilayer locations, predicted using the C16/DAcPC system, are in excellent agreement with available experimental data, in contrast to other surrogate systems.

Structural determinants of drug partitioning in n-hexadecane/water system.

Surrogate phases have been widely used as correlates for modeling transport and partitioning of drugs in biological systems, taking advantage of chemical similarity between the surrogate and the phospholipid bilayer as the elementary unit of biological phases, which is responsible for most of the transport and partitioning. Solvation in strata of the phospholipid bilayer is an important drug characteristic because it affects the rates of absorption and distribution, as well as the interactions with the membrane proteins having the binding sites located inside the bilayer. The bilayer core can be emulated by n-hexadecane (C16), and the headgroup stratum is often considered a hydrophilic phase because of the high water content. Therefore, we tested the hypothesis that the C16/water partition coefficients (P) can predict the bilayer locations of drugs and other small molecules better than other surrogate systems. Altogether 514 PC16/W values for nonionizable (458) and completely ionized (56) compounds were collected from the literature or measured, when necessary. With the intent to create a fragment-based prediction system, the PC16/W values were factorized into the fragment solvation parameters (f) and correction factors based on the ClogP fragmentation scheme. A script for the PC16/W prediction using the ClogP output is provided. To further expand the prediction system and reveal solvation differences, the fC16/W values were correlated with their more widely available counterparts for the 1-octanol/water system (O/W) using solvatochromic parameters. The analysis for 50 compounds with known bilayer location shows that the available and predicted PC16/W and PO/W values alone or the PC16/O values representing their ratio do not satisfactorily predict the preference for drug accumulation in bilayer strata. These observations indicate that the headgroups stratum, albeit well hydrated, does not have solvation characteristics similar to water and is also poorly described by the O/W partition characteristics.

Binding of Matrix Metalloproteinase Inhibitors to Extracellular Matrix: 3D-QSAR Analysis

Binding to the extracellular matrix, one of the most abundant human protein complexes, significantly affects drug disposition. Specifically, the interactions with extracellular matrix determine the free concentrations of small molecules acting in tissues, including signaling peptides, inhibitors of tissue remodeling enzymes such as matrix metalloproteinases, and other drug candidates. The nature of extracellular matrix binding was elucidated for 63 matrix metalloproteinase inhibitors, for which the association constants to an extracellular matrix mimic were reported here. The data did not correlate with lipophilicity as a common determinant of structure‐nonspecific, orientation‐averaged binding. A hypothetical structure of the binding site of the solidified extracellular matrix surrogate was analyzed using the Comparative Molecular Field Analysis, which needed to be applied in our multi‐mode variant. This fact indicates that the compounds bind to extracellular matrix in multiple modes, which cannot be considered as completely orientation‐averaged and exhibit structural dependence. The novel comparative molecular field analysis models, exhibiting satisfactory descriptive and predictive abilities, are suitable for prediction of the extracellular matrix binding for the untested chemicals, which are within applicability domains. The results contribute to a better prediction of the pharmacokinetic parameters such as the distribution volume and the tissue‐blood partition coefficients, in addition to a more imminent benefit for the development of more effective matrix metalloproteinase inhibitors.

Structure-based prediction of drug distribution across the headgroup and core strata of a phospholipid bilayer using surrogate phases.

Solvation of drugs in the core (C) and headgroup (H) strata of phospholipid bilayers affects their physiological transport rates and accumulation. These characteristics, especially a complete drug distribution profile across the bilayer strata, are tedious to obtain experimentally, to the point that even simplified preferred locations are only available for a few dozen compounds. Recently, we showed that the partition coefficient (P) values in the system of hydrated diacetyl phosphatidylcholine (DAcPC) and n-hexadecane (C16), as surrogates of the H- and C-strata of the bilayer composed of the most abundant mammalian phospholipid, PC, agree well with the preferred bilayer location of compounds. High P values are typical for lipophiles accumulating in the core, and low P values are characteristic of cephalophiles preferring the headgroups. This simple pattern does not hold for most compounds, which usually have more even distribution and may also accumulate at the H/C interface. To model complete distribution, the correlates of solvation energies are needed for each drug state in the bilayer: (1) for the H-stratum it is the DAcPC/W P value, calculated as the ratio of the C16/W and C16/DAcPC (W for water) P values; (2) for the C-stratum, the C16/W P value; (3) for the H/C interface, the P values for all plausible molecular poses are characterized using the fragment DAcPC/W and C16/W solvation parameters for the parts of the molecule embedded in the H- and C-strata, respectively. The correlates, each scaled by two Collander coefficients, were used in a nonlinear, mass-balance based model of intrabilayer distribution, which was applied to the easily measurable overall P values of compounds in the DMPC (M = myristoyl) bilayers and monolayers as the dependent variables. The calibrated model for 107 neutral compounds explains 94% of experimental variance, achieves similar cross-validation levels, and agrees well with the nontrivial, experimentally determined bilayer locations for 27 compounds. The resulting structure-based prediction system for intrabilayer distribution will facilitate more realistic modeling of passive transport and drug interactions with those integral membrane proteins, which have the binding sites located in the bilayer, such as some enzymes, influx and efflux transporters, and receptors. If only overall bilayer accumulation is of interest, the 1-octanol/W P values suffice to model the studied set.