METHOD OF MODELING AND PREDICTING BINDING OF LIGAND MOLECULES TO TARGET MOLECULES BY METHODS OF QUANTUM MECHANICS TAKING EFFECT OF SOLVENT INTO ACCOUNT
DESCRIPTION The invention relates to the development of new compounds for medicine, biotechnology, physics of condensed mediums and different fields of the science of materials, in which intermolecular interaction is significant, and may be used during the development of new medicaments by the in silico methods - i.e., by computational modeling a complex of a medicinal agent with a receptor or specific protein. Energetic values characterizing the intensity of the interaction of the protein and the medicinal agent may be obtained in the course of such modeling. A conclusion may be made on the basis of these values that relates to the stability of the complex of the medicinal agent with the receptor or specific protein. The stability of the complex is directly related to the possible effectiveness of the medicinal compound. Models of a force field are widely used at present for this purpose (for example, such as AMBER2002, OPLS-AA, CHARMM, ab-initio force field MMFF), in which the atoms in the molecules are considered as classical particles interacting with each other with the aid of interaction potentials. In spite of the wide use of force field models to obtain a numerical modeling of the interaction of molecules with each other, this methodology quite often inadequately describes the behavior of molecular systems. This is due to the fact that the behavior of atoms in molecules and electrons in atoms is subordinate to quantum laws and not to classical mechanics. At present quantum mechanics is widely used to calculate the characteristics of molecules of relatively small size (tens of atoms) and their interaction with each other. Recently, in view of the rapid increase of the power of computers, quantum calculations of molecular systems comprising hundreds and even thousands of atoms have become possible. The accuracy of and time for such calculations depend on the used methods of calculation and computer power. Taking the intensive development of the first and the rapid growth of the second into account, it is clear that the use of quantum mechanics for the numerical modeling of molecular systems will increase. In particular, in order to calculate the intensity of the binding of ligand molecules to a target molecule, the classical force field methods will be replaced by methods based on equations of quantum mechanics. However methods of quantum mechanics were earlier very rarely used for the description of biological molecules which, as a rule, are macromolecules comprising hundreds and thousands of atoms. Thus, the ab initio Hartree-Fock method and the density functional method (DFT) are used in papers [1-3] to describe the quantum part of a system which is
immersed in a molecular surrounding described classically with the aid of some force field. This is the so-called QM7MM approach in which the description of a part of the system takes place in accordance with the laws of quantum mechanics (QM), the other part of the system - in accordance with the laws of classical molecular mechanics (MM). In particular, the calculation of the energy of interaction (more exactly the enthalpy of interaction) between two biological molecules is carried out in the paper [2], in which an enzyme-substrate reaction is modeled. A simpler and accordingly more rapid quantum-chemical semiempirical method PM3 is used in the paper [4] to model the active center in a deaminase cytidine enzyme. The molecular system during the modeling contained an active center of 1330 atoms. In order to calculate such a large molecular system, the so-called Divide-and-Conquer (or abbreviated D&C or DAC) method was used. The semiempirical quantum-mechanical method AMI was used to study the electronic structure of biological macromolecules in a solution [5], wherein a continual model of a solvent was used. The enthalpies of the protein - protein and protein - DNA binding are obtained as the difference between the enthalpies of the formation of corresponding complexes and the enthalpies of the formation of separate components thereof, but the obtained data turned out to be sufficiently unrealistic (about 5-20 eV). Several organic macromolecules comprising from 256 to 9378 atoms were calculated within the frame of the semiempirical method PM3 with use of the D&C method within the frame of the package of the quantum-chemical program MOPAC [6]. Wherein, the calculations were carried out both in a vacuum and in an aqueous solution with the use of a COSMO approach, in which the solvent was continually modeled. Wherein, the energetic characteristics of the chemical systems being studied were not determined either, and the main attention was given to the attainment of optimized structures of the macromolecules and the time required for optimization. The majority of known quantum-mechanical studies of biological macromolecules are dedicated to calculation of the energetic characteristics of different enzymatic reactions. Wherein, the mechanisms of different reactions with the participation of water molecules, their transition complexes and corresponding energetic barriers were studied. However, in none of these works is there assessment of the intensity of the binding of proteins and ligands in their intermolecular complexes and the effect of the solvent on the processes under consideration in the macromolecules is not taken into account. Thus, it is evident from known literature that in application of quantum mechanics methods to modeling of biological macromolecules, there is no a single approach to the determination of the intensity of the binding of ligand molecules to target molecules. The structural and electronic characteristics of individual macromolecules with the use of
mechanisms and energetic barriers of the enzymatic reactions of disruption of the peptide bonds in proteins are mainly determined in the known methods with the aid of quantum mechanical methods. A comparison of such calculated data with experimentally obtained data shows unsatisfactory accuracy of the calculated values. Methods are not known from prior art for modeling energetic characteristics of the intermolecular interaction of ligand molecules and target molecules on the basis of quantum-mechanical calculations. Furthermore, in existing methods for quantum-mechanical studies of macromolecules, the effect of the solvent, which is very important during modeling of the binding of biological target molecules to ligand molecules, where hydrogen bonds play a significant role, is not taken into account in an appropriate manner. In view of this, the object of the instant invention is to create a method of modeling intermolecular interaction in an aqueous solution on the base of methods of quantum calculations, taking the effect of the solvent into account, both on a molecular level and in the form of a screening continuum, taking the spatial relaxation of atoms and the redistribution of electronic density during the formation of intermolecular complex into account. Such a method of modeling should be applied in the determination of the intensity of binding ligand molecules to target molecules independent of the makeupand size of these molecules, and also independent of the type of quantum calculation method used. This method should use the minimum number of adjusting parameters and be based as far as possible only on the main principles of quantum mechanics. The technical result is the creation of a method of modeling and predicting the binding of ligand molecules to target molecules by quantum mechanical methods taking the effect of the solvent into account. The indicated result is achieved in that a method of modeling and predicting the binding of target molecules and ligand molecules is proposed in which on the basis of the structure and makeup of the intermolecular target-ligand complex, the target molecule and the ligand molecule, models of an intermolecular complex and components thereof are formed, taking the solvent into account in an explicit and implicit manner, which models are further used during the modeling of the binding of the ligand molecule and the target molecule, taking the effect of the solvent into account by quantum mechanical methods with selected modeling parameters. Using the results of this modeling, the enthalpy, entropy and free energy of binding a selected ligand molecule and a target molecule, characterizing the intensity of binding the ligand molecule to the target molecule, are calculated, and for known ligands are compared with known experimental data. For new ligands the calculated values with selected parameters of modeling are a criterion of the intensity of binding the ligand molecule to the target molecule, and, consequently, are a criterion during the prediction of the ligands as potential medicinal substances. In the claimed method of modeling, the structure and makeup of the target molecule and
the ligand molecule are obtained by means of experimental measurements from the PDB (Protein Data Bank) database [7]. Furthermore, the structure, makeup and coordinates of the atoms of the new ligand molecule may be obtained with the aid of special programs - builders of new ligands. In the case where data on the Cartesian coordinates of hydrogen atoms are not present in the initial structures of the target molecules and ligand molecules, the coordinates of the missing atoms are added to the coordinates of the atoms of the target molecule and ligand molecule, using different computer programs for building molecular structures. The positions of the added hydrogen atoms in space are determined with molecular-mechanical or quantum- mechanical optimization under the condition at which the total energy of the target molecule and the ligand molecule reach the minimum during variation of the coordinates of the added hydrogen atoms. The structure of the intermolecular complex of a target molecule and ligand molecule is obtained with the aid of joining the ligand molecule and the target molecule or, as more often called, docking. In order to do this the coordinates of the atoms of these molecules are converted in such a manner that the least distance between each of the ligand atoms and predetermined atoms of the target molecule would be in a predetermined range. In order to take into account the effect of the solvent during modeling on the binding of the target molecule to the ligand molecule within the frame of the explicit model of the solvent, a set of coordinates of the solvent molecules are added to a set of coordinates of the intermolecular target-ligand complex and to a set of coordinates of the ligand molecule and the target molecule separately. This is done so that the solvent molecules themselves would fill, without mutual intersection or intersection with the atoms of the intermolecular complex, the target molecule and ligand molecule, all the space around these molecular systems. The position of the solvent molecules is optimized, as a rule, within the frame of the packets of molecular mechanics or molecular dynamics near the intermolecular complex and its components, taken individually. The quantum-mechanical models of the target molecule, ligand molecule and their intermolecular complex are built on the base of the obtained structures of that intermolecular complex and its components with a solvent taken into account explicitly. In order to do this the atoms of the intermolecular complex, target molecule, and ligand molecule together with solvent atoms are sorted into two groups. As regards the intermolecular complex, the first group should include all the atoms of the ligand molecule, and also atoms of the target molecule and the solvent molecule, where the distance from them to each of the atoms of the ligand molecule does not exceed a predetermined value, and the second group includes all the remaining atoms of the target molecule and all the remaining molecules of the solvent. As regards the model of the target molecule, the first group of atoms should include the same atoms of the target molecule which were included in the first group of atoms of the intermolecular complex. Furthermore, the
first group should include solvent molecules which are at a distance from the atoms of the target molecule included in the first group that does not exceed a predetermined value, and the second group includes all the remaining atoms of the target molecule and all the remaining solvent molecules. As regards the model of the ligand molecule, the first group of atoms should include all the atoms of the ligand molecule and also the solvent molecules which are at a distance from the atoms of the ligand molecule that does not exceed a predetermined value, and the second group includes all the remaining solvent molecules. As a rule, the atoms of the first group of all the molecular systems under consideration are modeled by methods of quantum mechanics, while atoms of the second group are modeled by methods of classical mechanics or atoms of the second group are withdrawn from the modeling process. In order to model the joint of two groups of atoms, calculated within the frame of quantum and classical methods, additional atoms are introduced near the region of the space where the atoms of the first and second groups are at a distance from each other that does not exceed a predetermined value. Wherein additional boundary atoms are calculated simultaneously in the classical and quantum approximations. In the case where only atoms of the first group are calculated within the frame of the quantum method, while atoms of the second group are withdrawn from the modeling process, the added atoms are also calculated within the frame of the quantum approach. Wherein, the additional atoms are added in such a manner that they form new covalent bonds in place of the disrupted covalent bonds between the atoms of the first and second groups. The charge states of different molecular groups are determined in the obtained models and in accordance therewith additional atoms are added to or removed from corresponding places in those molecular groups. The full charge of the built models of the intermolecular target-ligand complex, target molecule and ligand molecule is determined in accordance with those charge states. In the case where the solvent is water, then in order to better reproduce the formation and disruption of hydrogen bonds during the transition from the complex to its components, additional water molecules are introduced into the model. Molecules of water are added to the intermolecular complex so that they are located between the target molecule and the ligand molecule and form hydrogen bonds therewith. Molecules of water are added to the target molecule and to the ligand molecule so that with their participation the disrupted hydrogen bonds, which are present in the intermolecular complex between the target molecule and the ligand molecule and which were disrupted during the removal of the ligand molecule from the target molecule, are restored. However, in order to correctly take the effect of the solvent into account, it is not, as a rule, sufficient to explicitly take into account some amount of the solvent molecules, and in that
case it is necessary to take the solvent into account in an implicit model, where its screening effect is taken into account by introducing a continuous medium with the necessary value of dielectric permeability into the model. In the proposed method of modeling the bond between the target and the ligand, three types of models of an implicit account of the solvent are used, in particular, a model realizing solution of Poisson's equation, a model of a continuous conductor, and a model realizing solution of the Poisson-Boltzmann equation. All the data obtained at the step of building models of the target molecule, ligand molecule and intermolecular complex are used to form input files for quantum-mechanical modeling. These input files contain information on the coordinates and types of all the atoms of the molecular systems under consideration, their net charge and multiplicity, coordinates and type of solvent molecule atoms added in the explicit model to the intermolecular complex and its components, the parameters of the solvent in the implicit model. Furthermore, there is information in the input files on the parameters of the optimization process, in particular the method and parameters of the modeling method, type of optimization procedure, accuracy of self-consistency and minimization. The input files are used, so that with the aid of a corresponding program for modeling the total energy of a molecular system, based either on a combination of quantum-mechanical and classical approximation or only on quantum-mechanical approximation, to find the minimum of the total energy of the three molecular systems under consideration in the case where the positions of all or some of the additionally determined atoms are varied. The obtained values of the total energies of the target molecule, ligand molecule and their intermolecular complex are used to calculate the enthalpy of binding the components in the complex. The enthalpy of binding is calculated as the difference between the total energy of the intermolecular target-ligand complex in the solvent, on the one hand, and the sum of the total energies of the target molecule in the solvent and the ligand molecule in the solvent on the other hand while keeping or while not keeping a balance of the solvent molecules. In accordance with the instant invention, calculation of the enthalpy of binding the target molecule and the ligand molecule may be carried out using another method in which the total energies of the target molecule, the ligand molecule and their intermolecular complex are calculated not in the minimum point, but as the average values of the total energies obtained for a number of configurations of the intermolecular target-ligand complex. Wherein the coordinates of the atoms of the target molecule, ligand molecule and their intermolecular complex, corresponding to each of such configurations, may be obtained within the frame of molecular- dynamic modeling. In such a type of modeling the structure of the target molecule and the ligand molecule and their intermolecular complex, obtained by joining the coordinates of its components, is also used as the initial data. In the case of molecular-dynamic modeling, both
explicit and implicit models of the solvent, which are introduced by the methods described earlier, may be used. In the case of an intermolecular complex, independent of the model of the solvent, modeling is carried out by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the target molecule and the ligand molecule, and taking into account the presence of the solvent molecules, wherein the state of the molecular system at which it reaches thermal balance is determined. Then the trajectory of movement of all the atoms of the system after they have reached thermal balance is determined and the coordinates of all the atoms of the interacting target molecule, ligand molecule and solvent molecules are determined after predetermined intervals of time a certain number of times for a predetermined set of configurations of the intermolecular target-ligand complex. Similarly, a set of configurations of the target molecule and the ligand molecule is obtained. For each set of coordinates of the intermolecular target-ligand complex, target molecule, ligand molecule in the explicit or implicit model of the solvent, all the atoms are sorted into groups in a manner similar to that carried out earlier. An input file for quantum-mechanical modeling of the intermolecular target-ligand complex in the solvent, is formed, as described above, for each configuration of the intermolecular target-ligand complex, target molecule, ligand molecule. For each input file, i.e., for each configuration of the intermolecular target- ligand complex, target molecule, ligand molecule, calculation of the total energy of these molecular systems in the solvent is carried out with the aid of the selected method of modeling and without optimization of the spatial structure of the complex. Then the average values of the total energy of the intermolecular target-ligand complex, target molecule, ligand molecule in the solvent are calculated according to all the selected set of configurations of the molecular system. Using the obtained average values of the total energies of the intermolecular target-ligand complex in the solvent, the target molecule in the solvent and the ligand molecule in the solvent, the enthalpy of the binding of the ligand molecule to the target molecule in the solvent is calculated as the difference between the total energy of the intermolecular target-ligand complex in the solvent, on the one hand, and the sum of the total energies of the target molecule in the solvent and the ligand molecule in the solvent, on the other hand, while keeping a balance of solvent molecules. Calculation of the entropy of binding the ligand molecule to the target molecule in the binding solvent may be divided into calculation of three components thereof, in particular, vibrational, translational and rotational components. The entropy, related to a loss of the vibration degrees of freedom during the transition from free unbound states of the target molecule and the ligand molecule in the solvent to their intermolecular complex, is determined by the vibration frequencies of the molecules
participating in the process under consideration. In order to obtain the oscillating frequencies of the molecular systems under consideration, the force constants or the second derivatives of the total energies, relating to the intermolecular complex, the target molecule and ligand molecule, are calculated within the frame of the used method of quantum-mechanical modeling in respect to all the coordinates of the atoms in the points of corresponding rninimums. A contribution to the change of entropy during binding of the ligand molecule to the target molecule is made by the fact that these molecules during binding lose both translational and rotational degrees of freedom, which they had in the free unbound state. Only the coordinates and atomic weights of the atoms of the target molecule, ligand molecule and their intermolecular complex are used for calculation of this part of the entropy, and this calculation, as a rule, is carried out outside the frame of the quantum-mechanical modeling. The complete change of the entropy when the ligand molecule is bound to the target molecule is calculated as the sum of changes of different components thereof. The process of binding the ligand molecule to the target molecule takes place in a solvent. Each of the molecules creates a cavity in the solvent, pushing aside the solvent molecules. This process requires the expenditure of a certain free energy, which is called cavitation energy, and it depends only on the properties of the solvent and the form of the molecule disposed in the solvent. The change of the cavitation energy during binding of the target molecule and the ligand molecule is additionally calculated and is taken into account during the calculation of the free binding energy in the complex. The latter is calculated as the sum of the enthalpy, entropy members and the member characterizing a change of the cavitation energy. The calculated enthalpy, entropy and free energy of binding are energetic characteristics of the intensity of binding the target molecule and the ligand molecule and are used to predict the binding of new ligands to target molecules. The invention is explained by examples of embodiment, illustrated by drawings in which the following is presented: Fig. 1 shows a flowchart of a method of quantum-mechanical (QM) modeling and predicting the binding of a ligand molecule (L) to a target molecule (T); Fig. 2 shows a flowchart of a method of selecting a target molecule (T) and a ligand molecule (L) and building an intermolecular complex consisting of the ligand and target (L-T); Fig. 3 shows a flowchart of a method of building quantum-mechanical (QM) models of an intermolecular ligand-target complex (L-T) and components thereof with explicit and implicit accounting of the solvent; Fig. 4 shows a flowchart of a method of quantum-mechanical (QM) modeling of intermolecular interaction between a ligand molecule (L) and a target molecule (T), taking the
solvent into account, including calculation of the enthalpy, entropy, free energy of binding the ligand (L) and target (T); Fig. 5 shows a flowchart of a molecular-dynamic (MD) method of quantum-mechanical (QM) modeling of binding ligand molecules (L) to target molecules (T) in a solvent; Fig. 6 shows the structural formulas (left column) and the quantum-mechanically optimized structures (right column) of ligands corresponding to ligand-protein complexes, which were selected as test systems for modeling binding in a ligand-protein system by methods of quantum mechanics; Fig. 7 shows a model of an active center of a 1AJ6 complex (system 5 according to Table 1). The protein is presented by thin sticks, the water molecules are presented by thick sticks, the ligand molecule is presented by balls and sticks. The dashed lines designate hydrogen bonds in a ligand-protein complex; Fig. 8 shows a scheme of a method WI of adding water molecules to components of a ligand-protein complex; Fig. 9 shows a scheme of a method W2 of adding water molecules to components of the ligand-protein complex. Fig. 1 shows a general scheme of the claimed method of modeling and predicting the binding of a ligand molecule to a target molecule. In step 1 a target molecule and a ligand molecule are selected and their structure is obtained. In step 2 an intermolecular complex is formed from the selected target molecule and ligand molecule. On the basis of the structure and the makeup of that intermolecular complex, in step 3 models of the intermolecular target-ligand complex and its components are formed with explicit and implicit accounting of the solvent, which are further used in step 4 during the modeling of binding the ligand molecule and the target molecule, taking into account the effect of the solvent by methods of quantum mechanics with selected modeling parameters. Using the results of this modeling, numerical data are calculated in step 5 that characterize the intensity of binding the ligand molecule to the target molecule, in particular, the enthalpy, entropy and free energy of binding the selected ligand molecule and target molecule, and these values are recorded in the form of numerical data on an information carrier. In step 6, similar values characterizing the intensity of binding the predetermined ligand molecule to the predetermined target molecule are experimentally measured. Then the values calculated during the modeling and characterizing the intensity of binding the ligand molecule to the target molecule are compared in step 7 with corresponding values measured in step 6. In the case where the calculated and measured values differ from the predetermined limits, in step 8 the model of building the intermolecular complex and the parameters of the quantum-mechanical modeling of the interaction of the target molecule and the ligand molecule are corrected, in particular, the model of the account of the effect of the solvent
on the binding is corrected. Then modeling is carried out again in step 4 and the values characterizing the intensity of binding the ligand molecule to the target molecule are calculated in step 5. This process of developing the optimum model and parameters of modeling is continued until the difference between the calculated and measured values, characterizing the binding of the molecules to each other, becomes less than the predetermined value. After that the process of adjusting the parameters of modeling for the predetermined target molecule and ligand molecule pair is considered to be finished. Then, another ligand molecule and target molecule pair is selected and the steps described above are repeated for this pair, i.e., the steps of forming a model of the intermolecular complex (3), quantum-mechanical modeling (4), calculating values characterizing the binding of the ligand molecule to the target molecule (5), measuring these values (6), comparing the calculated and measured values (7) and correcting the parameters (8). The values of the modeling parameters, i.e., values that are obtained for the preceding target molecule and ligand molecule pair, are used as the initial parameters. Then the next target molecule and ligand molecule pair is selected and so on. As a result of such a self-consistent procedure, parameters of modeling are determined, which make it possible to obtain with predetermined accuracy meanings of values characterizing the intensity of binding a ligand molecule to a target molecule, in particular, enthalpy, entropy and free energy of binding, for a certain plurality of target molecule and ligand molecule pairs. This process of expanding the plurality of target molecule and ligand molecule pairs and self- consistent adjustment of the parameters of modeling continues until for the next newly-taken target molecule and ligand molecule pair the values, characterizing the intensity of binding and calculated with the use of modeling parameters obtained for other pairs of molecules, will not coincide with predetermined exactness with the values measured for this pair of molecules and characterizing the intensity of binding. After that the final values of the parameters of modeling are recorded in step 7 (Fig. 1) on an information carrier and further they are used for prediction. In step 9 affirmation of the correctness of the selected methods of forming the quantum- mechanical models and methods of modeling the intermolecular interaction in the complex is carried out by demonstrating coincidence with the predetermined accuracy (which should of the order of the accuracy of the quantum calculations 0.1-0.2 kcal/mol themselves) of the values experimentally obtained and calculated within the frame of the developed model and characterizing the intensity of binding for the next target molecule and ligand molecule pairs, which were not used in the iterative procedure of the development of the method of modeling. Further the procedure for building the quantum-mechanical model and modeling the intermolecular interaction in the complex consisting of a ligand molecule and a target molecule is used to determine in step 10 the group of new ligands that have better binding to a
predetermined target molecule. Predicting in step 11 the ligand molecules which sufficiently strongly bind to a predetermined target molecule and may be new potential medicinal substances consists of the following. A plurality of ligand molecules is formed for a predetermined target molecule, which are candidates for predicting, on the basis of different databases and libraries of fragments with the aid of programs using certain rules for building ligands in the active center of the target molecule. The initial databases comprise information on the structure of the molecules or the structure of the fragments of molecules, which information is obtained in different experiments, or comprises the coordinates of atoms of the molecules or fragments of the molecules, built with the aid of different computer programs. The number of molecules in such a plurality of ligand molecules, which are candidates for predicting, may be huge, up to several millions, and each of these molecules exists in reality or may be synthesized. It is important that information on the chemical makeup, three-dimensional structure and coordinates of all the atoms be provided for the molecules of that plurality. For each ligand molecule among the selected plurality, modeling of the process of binding is carried out in step 5 (Fig. 1) and in step 6 calculation of the values characterizing the binding, taking into account the effect of the solvent, the ligand molecule with a predetermined target molecule, in particular the enthalpy, entropy and free energy of binding with use of the final values of the modeling parameters obtained at the preceding step. As a result of such modeling, values are obtained that characterize the intensity of binding all the ligand molecules from the selected plurality to the predetermined target molecule. According to the results of modeling, a number of the most promising ligand molecules are selected that have the best meanings of the values characterizing the intensity of the binding. Those ligand molecules for which good binding is predicted may be used in real chemical experiments in respect to binding to a predetermined target molecule. It is obvious that with such an approach the number of ligand molecules in respect to which experimental testing should be carried out is significantly reduced. Thus, the proposed method of computer prediction of the ligand molecules which may be new potential medicinal substances makes it possible to avoid significant time and financial expenditures for synthesis of new substances and conduction of experiments on the measurement of values characterizing the intensity of binding ligand molecules to target molecules. Further, a comprehensive description is provided of the main steps, presented in Fig. 1, of the claimed quantum-mechanical method of modeling the binding of a ligand molecule to a target molecule, taking the effect of a solvent into account. Fig. 2 shows a flowchart illustrating a method of forming an intermolecular complex of a target molecule and a ligand molecule (step 2 in Fig. 1) in accordance with the invention. On the basis of the plurality 12 of known ligand molecules in step 13 a target molecule is
selected for which it is necessary to carry out modeling and predicting the binding thereto of different ligand molecules. The structure and makeup of this target molecule is obtained in step 14 with the aid of experimental measurements, for example, with the aid of X-ray structural analysis, by methods of scattering neutrons or of nuclear magnetic resonance. The structure and makeup of the target molecule may also be obtained in the form of a set of Cartesian coordinates of corresponding atoms obtained in the process of molecular-mechanical or molecular-dynamic modeling, carried out after reading from corresponding information carriers preliminarily stored data extracted from different databases on the structure of a target molecule. In the case where the target molecule is a protein, its structure and makeup are obtained from the database of the PDB (Protein Data Bank) [7], accessible via the Internet. As a rule, data on the Cartesian coordinates of hydrogen atoms is not present in the initial structures of target molecules. The coordinates of the missing hydrogen atoms are added in step 15, there where it is necessary, to the coordinates of heavy atoms of the target molecule, using different computer programs for building molecular structures. The optimized meanings of the Cartesian coordinates of the added hydrogen atoms are determined in step 16 from the condition in accordance with which the total energy of the target molecule reaches the minimum value (in the case of molecular-mechanical or quantum-mechanical optimization) while varying the coordinates of the added hydrogen atoms. Information on the structure and the makeup, and also on the coordinates of all the atoms, including the coordinates of the added hydrogen atoms of the target molecule, is recorded in step 17 on an information carrier in a corresponding format. On the basis of the plurality 18 of known ligand molecules in step 19 a ligand molecule is selected in respect to which modeling the binding of it to a selected target molecule will be carried out. Besides that, a plurality of 20 new ligands may be built on the basis of the different kind of databases or libraries of fragments of ligands, and a ligand molecule may be selected in step 19 from that plurality, in respect to which the prediction of its binding to a selected target molecule will be subsequently carried out. The structure, makeup and coordinates of the atoms of a known ligand molecule are obtained in step 21 in a manner similar to that done in respect to the target molecule, in particular by or with the aid of experimental measurements, within the frame of the molecular-mechanical or molecular-dynamic modeling or on the basis of the PDB database. The structure, makeup and coordinates of the atoms of the new ligand molecule may be obtained with the aid of special programs - builders of new ligands. In the case where the experimentally known or newly formed structure of the ligand molecule does not comprise hydrogen atoms, then in step 22, where necessary, coordinates of the missing hydrogen atoms are added to the coordinates of the heavy atoms of the ligand, using different computer programs for forming molecular structures. The coordinates of these
hydrogen atoms are determined in step 23 by minimizing (either within the frame of molecular mechanics or within the frame of quantum mechanics) the total energy of the ligand molecule while varying the coordinates of the added hydrogen atoms. The coordinates of all the atoms of the ligand molecule, including the coordinates of the added hydrogen atoms, are recorded in step 24 on an information carrier in a corresponding format. Then in step 25 by joining the ligand molecule and the target molecule a structure of an intermolecular complex of a target molecule and a ligand molecule (2 in Fig. 1) is obtained, which is subsequently used in quantum-mechanical modeling. In order to do this the coordinates of all the atoms of the target molecule and the ligand molecule are read, including the coordinates of the added hydrogen atoms, from the information carrier into the computer memory. After that a transformation of the coordinates of the atoms of these molecules is carried out in such a manner that the smallest distance between each of the ligand atoms and the predetermined atoms of the target molecule is within a predetermined range. Such a procedure of joining or. as it is more often called, docking, is necessary for the overlapping of the two molecules in space, since usually the coordinates of their atoms are obtained from different sources of information. The thus obtained coordinates of the atoms of the intermolecular complex of the target molecule and the ligand molecule are recorded in step 26 into the computer memory and used to build a quantum-mechanical model of a corresponding intermolecular complex. Further, the intermolecular complex and its components are converted in such a manner that the effect of the solvent on the binding of the target molecule and the ligand molecule is taken into account on the basis of a model of explicit account of the solvent. In order to do this in step 27 a plurality of coordinates of solvent molecules are added to a plurality of coordinates of the intermolecular target-ligand complex. In steps 28 and 29 the plurality of coordinates of solvent molecules are added to the plurality of coordinates of the components of this complex, taken individually, in particular to the plurality of coordinates of the target molecule (in step 28) and to the plurality of coordinates of the ligand molecule (in step 29). In the case of the intermolecular complex this is done in step 27 by computer calculation of the coordinates corresponding to the positions in space of the plurality of solvent molecules, so that their average density would correspond to the experimentally observed density of the solvent. The solvent molecules should fill, without mutually intersecting or intersecting with atoms of the intermolecular complex of the target molecule and the ligand molecule, the space both between the target molecule and the ligand molecule and around them. Similarly, the coordinates of the plurality of solvent molecules are added only to the plurality of coordinates of the target molecule in step 28, so that their average density would correspond to the experimentally observed density of the solvent, and they would fill without mutual intersection or intersection
with atoms of the target molecule the space around it. Also, the coordinates of the plurality of solvent molecules are added only to the plurality of coordinates of the ligand molecule in step 29 so that their average density would correspond to the experimentally observed density of the solvent, and they fill without mutual intersection or intersection with the atoms of the ligand molecule the space around it. The position of the solvent molecules is optimized, as a rule, within the frame of the packages of molecular mechanics or molecular dynamics near the intermolecular complex in step 30, and also near its components taken individually, in particular, near the target molecule in step 31 and the ligand molecule in step 32. In step 33 the coordinates of the atoms of the intermolecular target-ligand complex: together with the coordinates of the atoms of the solvent molecules near the surface of the complex are recorded on the information carrier in a corresponding format. Furthermore, in step 33, the coordinates of atoms of the target molecule together with the coordinates of the atoms of the solvent molecules, located near its surface, and the coordinates of the atoms of the ligand molecule together with the coordinates of the atoms of the solvent molecules, near its surface, are recorded individually on the carrier in a corresponding format. Then, on the basis of the structure of the obtained intermolecular complex, consisting of a target molecule and a ligand molecule, and also molecules of the solvent, quantum-mechanical models of that intermolecular complex and its components are formed (step 3 in Fig. 1) with explicit and implicit account of the solvent and used in quantum-mechanical modeling. Fig. 3 shows a flowchart illustrating the method according to the invention of forming such quantum- mechanical models and, accordingly, forming input files that are subsequently used for modeling the binding in the intermolecular target-ligand complex. In step 34 with the aid of a computer program the coordinates of the atoms of the intermolecular target-ligand complex together with the added solvent molecules are read from the computer memory, as are the coordinates of the target molecule atoms together with the solvent molecules added thereto and the coordinates of the atoms of the ligand molecule together with the solvent molecules added thereto. Then the distances from the atoms of the ligand molecule in the intermolecular complex to each atom of the target molecule and the solvent are calculated. In order to form a model of the intermolecular target-ligand complex, in step 35 the atoms of this complex are sorted into two groups. This is done in such a manner that the first group includes all the atoms of the ligand molecule and also the atoms of the target molecule and the solvent molecule, the distances from which to each of the atoms of the ligand molecule do not exceed a predetermined value, and the second group includes all the remaining atoms of the target molecule and all the remaining molecules of the solvent. Wherein, if even one atom of any solvent molecule is within the limits of the predetermined distance from the atoms of the ligand molecule, then all the atoms of that solvent molecule are included in the first group of atoms.
After that in that same step 35 an input file is begun to be formed with use of the coordinates of the atoms of the first and second groups, the input file being used for calculation of the intermolecular target-ligand complex in a solvent with the aid of a corresponding program for modeling in such a manner that the atoms of the first group are modeled by methods of quantum mechanics, and the atoms of the second group are modeled by methods of classical mechanics. Wherein the operator of the total energy is built in the form of the sum of the following operators: an operator of the total energy for the quantum part, consisting of atoms of the first group, an operator of the total energy for the classical part, consisting of atoms of the second group, and an operator of the interaction between the classical and the quantum parts. In order to form a model of the target molecule, in step 36 the atoms of the target molecule and relating thereto molecules of the solvent are sorted in such a manner that a first group of atoms includes the same atoms of the target molecule which are included in the first group of atoms that is formed in step 35. Furthermore, included in the first group are also solvent molecules which are at a distance from the atoms of the target molecule that are included in the first group that does not exceed a predetermined value, and in the second group are included all the remaining atoms of the target molecule and all the remaining solvent molecules. After that in that same step 36 an input file is begun to be formed with use of the coordinates of the atoms of the first and second groups, the input file being used for calculation of the target molecule in the solvent with the aid of a corresponding program for modeling in such a manner that the atoms of the first group are modeled by methods of quantum mechanics, and the atoms of the second group are modeled by methods of classical mechanics. In order to form a model of the ligand molecule, in step 37 the atoms of the ligand molecule and the molecules of the solvent are sorted in such a manner that a first group of atoms includes all the atoms of the ligand molecule and also the solvent molecules which are at a distance from the atoms of the ligand molecule that does not exceed a predetermined value, and in the second group are included all the remaining solvent molecules. After that in that same step 37 an input file is begun to be formed with use of the coordinates of the atoms of the first and second groups, the input file being used for calculation of the ligand molecule in the solvent with the aid of a corresponding program for modeling in such a manner that the atoms of the first group are modeled by methods of quantum mechanics, and the atoms of the second group are modeled by methods of classical mechanics. In the case where the target molecule is a protein, sorting the atoms into two groups (in step 35) in the intermolecular complex is carried out in the following manner. The first group of atoms includes atoms of the ligand molecule, atoms of the solvent molecules, which are at a distance from the atoms of the ligand molecule that does not exceed a predetermined value, and atoms of the target molecule, which form whole amino acid residues, at least one atom of which
is at a distance from each atom of the ligand molecule that does not exceed a predetermined value. The second group includes all the remaining atoms of the target molecule and the solvent molecule. During the building of a model of the target molecule, which is a protein, included in the first group of atoms (in step 36) are atoms of the target molecule, which compose whole amino acid residues, relating to the first group of atoms for the intermolecular complex. For more correct modeling of the joining of the two groups of atoms for which processing of the data is carried out within the frames of the quantum and classical methods, in steps 38, 39 and 40 the coordinates of additional atoms are added into the input files for a corresponding program of modeling the intermolecular ligand-target complex, the target molecule and the ligand molecule respectively. These additional atoms are positioned near the region of the space where the atoms of the first and second group are located at a distance from each other that does not exceed a predetermined value. Modeling the binding of components in the intermolecular target-ligand complex may also be carried out in another manner so that the coordinates of the atoms of the second group are not included in the input files of a corresponding program of modeling the intermolecular target-ligand complex, the target molecule and the ligand molecule, and the atoms of the first group are modeled by methods of quantum mechanics. Wherein, the coordinates of the additionally added atoms are introduced into the input files of a corresponding modeling program so that they together with the atoms of the first group are modeled by methods of quantum mechanics. With such a method of modeling, in the case where unsaturated broken bonds remain in atoms of the first group, which had covalent bonds with atoms of the second group, it is necessary to restore or compensate loss of these covalent bonds. In that case, during the formation of a model of an intermolecular target-ligand complex in step 41 and a model of a target molecule in step 42, additional atoms are added to the atoms of the first group in such a manner that they form new covalent bonds in place of the broken covalent bonds between atoms of the first and second groups. In the case where the target molecule is a protein, the broken covalent bonds of the N- ends of the main chain of the protein are closed by -C(O)-(CH
3) groups, while the broken covalent bonds of the C-ends of the main chain of the protein are closed by — N(H)-(CH
3) groups. During the formation of input files for modeling the binding in the ligand - protein complex in steps 43, 44 and 45, the charge states of different molecular groups of the intermolecular target - ligand complex, target molecule and ligand molecule in the solvent are determined (as a rule, this is determined by the physicochemical properties of these molecular groups in the solvent or directly in experiments). In accordance with this, additional atoms are added or withdrawn in corresponding places of the intermolecular target - ligand complex in step 46, the target molecule in step 47 and the ligand molecule in step 48. When the solvent is
water, the required charge states of different molecular groups of components of the intermolecular target-ligand complex, target molecule and ligand molecule in the solvent are obtained by the addition or withdrawal of additional hydrogen atoms in corresponding places of the target molecule and the ligand molecule. Further, in step 49 determination of the full charge of the intermolecular target-ligand complex, target molecule and ligand molecule is carried out. The full charge of the intermolecular complex and its components separately, i.e., target molecule and ligand molecule, is obtained by a simple addition of all the charges which were determined for separate molecular groups included in the makeup of the molecular systems under consideration. The obtained data on the full charges of the models of the intermolecular target-ligand complex, target molecule and ligand molecule are further used during the formation of three input files for a corresponding program of modeling binding in the complex together with the multiplet data on the molecular systems under consideration, the coordinates and types of their atoms. In the case where the solvent is water, in order to more exactly take into account in an explicit form the effect of the solvent on the process of binding the target molecule and the ligand molecule, it is necessary to add the coordinates of the atoms of additional water molecules into each of the three input files of the target molecule, ligand molecule and their intermolecular complex. These additional water molecules are added in the process of modeling in order to better reproduce the formation and disruption of hydrogen bonds, which play a sign ificant role in the intermolecular target - ligand complex of the target molecule and ligand molecule. In step 50, additional water molecules are added to the intermolecular complex so that they would be between the target molecule and the ligand molecule and form therewith hydrogen bonds. The addition of water molecules to the target molecule in step 51 and to the ligand molecule in step 52 is carried out so that with their participation the disrupted hydrogen bonds, which were present in the intermolecular complex between the target molecule and the ligand molecule and which were disrupted during the withdrawal of the ligand molecule from the target molecule, are restored. Furthermore, in the case where the water molecules are arranged between the target molecule and the ligand molecule and were bonded to them both by hydrogen bonds, these water molecules are added twice each and to the input file for calculation of the target molecule in step 5 1 and to the input file for calculation of the ligand molecule in step 52 during the separate modeling of the components of the intermolecular complex. As a rule, the solvent has a value of dielectric permeability that is different from the value of dielectric permeability of a vacuum and the value of dielectric permeability of a target macromolecule. For example, the dielectric permeability of water is equal to approximately 80 at room temperature and normal pressure, while the dielectric permeability of a protein dissolved in
water does not exceed several units, and wherewith the dielectric permeability of a vacuum is equal to 1. In that case the presence of a solvent around the molecules under consideration significantly changes the electrical potential and the field, which are created by the charges of the atoms around the molecules. In order to take this screening effect of the solvent into account, it is not, as a rule, sufficient to clearly take an explicit amount of the solvent molecules into account, and in that case it is necessary to take the solvent into account in an implicit model, when its screening effect is taken into account by introducing a continuous medium with the necessary value of dielectric permeability into the model. Within the frame of the implicit model of the solvent, in step 53 the parameters of the solvent in an implicit model are introduced into the input files for calculation of the intermolecular complex, the target molecule and ligand molecule with the aid of a corresponding program of modeling molecular systems in a solvent. In order to do this the coordinates of points on the surface accessible by the solvent are formed with use of the coordinates of the atoms of the target molecule, ligand molecule and solvent molecules, the coordinates of which are introduced into corresponding input files,. This surface encompasses the target molecule, ligand molecule and solvent molecules in the case of calculation of the intermolecular complex in the solvent or encompasses only the target molecule and the solvent molecules in the case of calculation of the target molecule in the solvent or encompasses only the ligand molecule and the solvent molecules in the case of calculation of the ligand molecule in the solvent. A mathematical model of the implicit account of the solvent is formed with the aid of points that are on those surfaces. Three types of a model of the implicit account of the solvent are given consideration in the proposed invention. In the first of these a combination of software realizing the solution of a Poisson equation and describing the interaction of the target molecule and the ligand molecule with a continuous medium filling the whole space outside the surface accessible to the solvent and having a dielectric permeability equal to the dielectric permeability of the solvent, is used as a mathematical model of the implicit account of the solvent. In the case of use of a solvent with a large value of the dielectric permeab ility of the solvent (for example, water has a dielectric permeability of about 80 at room temperature), a combination of software realizing the solution of equations describing the interaction of a target molecule and a ligand molecule with a continuous conductor filling all the space outside the surface accessible to the solvent is used as a mathematical model of the implicit account of the solvent. In the case of use of a solvent characterized by a dielectric permeability and conductivity differing from zero, a combination of software realizing a solution of the Poisson-Boltzmann
equation is used as a mathematical model of the implicit account of the solvent. This equation describes the interaction of the target molecule and the ligand molecule with a continuous medium, having a dielectric permeability differing from zero and a conductivity differing from zero, and filling all the space outside the surface accessible to the solvent. The formation of the files describing the models for quantum-mechanical calculation of the intermolecular target-ligand complex, target molecule and ligand molecule in a solvent, ended at step 54, and the corresponding input files for the quantum-mechanical calculation are recorded on an information carrier. As a rule these files are used, depending on the selected method of calculation, in two variants of the quantum-mechanical modeling. In the first variant of modeling, calculation of the intermolecular complex, the target molecule, the ligand molecule in a solvent is carried out so that a first group of atoms of the complex, the target molecule and the ligand molecule (described above) is calculated by methods of quantum mechanics, and a second - by methods of classical mechanics. In the second variant of modeling, calculation of the intermolecular complex, the target molecule, the ligand molecule in a solvent is carried out so that the first group of atoms is calculated by methods of quantum mechanics and the second is not included in the modeling. Further, a description is provided of the process of quantum-mechanical modeling of the interaction of the ligand molecule and the target molecule, taking the solvent into account (step 4 in Fig. 1), and calculation of the energetic values characterizing the binding of the target molecule and the ligand molecule (step 5 in Fig. 1). In Fig. 4 a flowchart is presented that illustrates the method, conforming with the invention, of quantum-mechanical modeling of the intermolecular interaction between a ligand molecule and a target molecule, taking the solvent into account, including calculation of the enthalpy, entropy and free energy of binding the ligand and target. Input files in which information is recorded on the structure and properties of the intermolecular complex, the target molecule and the ligand molecule, are read from the information carrier in step 55. Those input files contain information on the coordinates and types of all the atoms of the molecular systems under consideration, their complete charge and multiplet properties, coordinates and type of atoms of the solvent molecules added in the explicit model to the intermolecular complex and its components, parameters of the solvent in the implicit model. Furthermore, the input files contain information on the parameters of the optimized process, in particular the method and parameters of quantum-mechanical modeling, the type of optimized procedure, the frame of self-agreement and minimization. In step 56, different methods of quantum-mechanical modeling of the total energy of the molecular systems may be selected. In the case of the approach when a portion of the atoms of
the intermolecular complex, the target molecule and the ligand molecule is calculated within the frame of quantum-mechanical approximation, and the other part of the atoms is calculated within the frame of classical approximation, different combinations of classical molecular-mechanical and molecular-dynamic methods are used on the one hand, and quantum-mechanical methods on the other hand. The latter may vary from strictly non-empiric methods and density functional methods [1-3] to more rapid, but less exact semiempirical methods [4-7]. Such an approach, becoming widely known under the designation QM/MM, is realized at the present time in the form of a number of applied program packets [8]. Furthermore, another approach is also possible, in which a part of the atoms of the intermolecular complex, target molecule and ligand molecule is calculated within the frame of quantum-mechanical approximation, while another part of the atoms is completely withdrawn from the modeling process. In order to more exactly reproduce the energetic characteristics of binding in the intermolecular complex, in the case of such an approach it is necessary to take as many as possible sections of target molecules, describing the place of binding to the ligand molecule with minimum influence of edge effects, for the quantum-mechanical modeling. Wherein the most promising methods of quantum-mechanical modeling are the relatively rapidly calculated semiempirical methods, which have recently been realized in a number of applied packets of programs with accelerated count [9]. However, the use of more exact, but slower non- empirical methods and methods of density functional, is also possible. At the beginning of the optimization procedure in steps 57, 58 and 59, determination of variable parameters of the intermolecular target-ligand complex, the target molecule, the ligand molecule and relating to those systems solvent molecules is carried out, in particular, determination of the variable coordinates of atoms, at which the minimum of the total energy of each of the molecular systems will be during a change thereof, is carried out Then in step 60, with the aid of a corresponding program for modeling the total energy of the molecular system, based either on a combination of quantum-mechanical and classical approximation or based only on quantum-mechanical approximation, the minimum is found of the total energy of the intermolecular target-ligand complex and the solvent with variation of the positions of all or some predetermined atoms of the complex and the solvent. After that the found new coordinates of the atoms and the total energy of the system, consisting of the intermolecular target-ligand complex and the solvent, and corresponding coordinates of atoms, are recorded on the information carrier. Similarly, in step 61, the minimum of the total energy of the target molecule and the solvent during variation of the positions of all or a part of the atoms of the target molecule and the solvent is found, and the found new coordinates of the atoms and total energy of the system, consisting of the target molecule and solvent, and of corresponding coordinates of the atoms are
recorded on the information carrier. Finally, in step 62, the minimum of the total energy of the ligand molecule and the solvent during variation of the positions of all or a part of the atoms of the ligand molecule and the solvent is found, and the found new coordinates of the atoms and total energy of the system, consisting of the ligand molecule and solvent, and of corresponding coordinates of the atoms are recorded on the information carrier. Using the obtained total energies of the intermolecular target - ligand complex and solvent, the target molecule and solvent, and also the ligand molecule and solvent, in step 63 the enthalpy of binding the ligand molecule to the target molecule in the solvent is calculated. Wherein, two cases of calculation of the enthalpy are possible. The first case is realized when a balance of the solvent molecules is maintained in the process of binding the ligand molecule to the target molecule, i.e., in the case where the sum of the solvent molecules included in the calculation of the target molecule and the solvent molecules included in the calculation of the ligand molecule is equal to the number of solvent molecules included in the calculation of the intermolecular target - ligand complex. In that case in step 63, the enthalpy of binding ΔH
bi
nd. the target molecule and the ligand molecule is calculated as the difference between the total energy of the intermolecular target - ligand complex in the solvent, on the one hand, and the sum of the total energies of the target molecule in the solvent and the ligand molecule in the solvent, on the other hand, in accordance with the equation: unbin .
_ ϋcomplex_NsolventmoI. - (, arget_Lsolveπtmcιl. ' ■'^ligand_MsoIventmol. wherein E
oompι
ex_Nsoive
nt
moi. is the total energy of the target-ligand complex, bound to N molecules of the solvent, E
targ
et Lsoi
ve
nt
moi. is the total energy of the target molecule bound to L molecules of the solvent, Eij
ga„
Msoi
veπ
tmoi. is the total energy of the ligand molecule bound to M molecules of the solvent. In the case where balance of the solvent molecules is observed N = L + M. The second case is realized when a balance of the solvent molecules is not maintained in the process of binding the ligand molecule to the target molecule, i.e., in the case -where the sum of the solvent molecules included in the calculation of the target molecule and the solvent molecules included in the calculation of the ligand molecule is not equal to the number of solvent molecules included in the calculation of the intermolecular target-ligand complex. In that case in step 63, the enthalpy of binding ΔH
bi
πd. the target molecule and the ligand molecule is also calculated as the difference between the total energy of the intermolecular target-ligand complex in the solvent, on the one hand, and the sum of the total energies of the target molecule in the solvent and the ligand molecule in the solvent, on the other hand. HoΛvever, during calculation of the enthalpy of binding ΔH
bi„d. in that step, the energy should be taken into
account which it is necessary to spend in order to remove from the solvent that amount of solvent molecules which it is necessary to add into the system under consideration in order to observe a balance in respect to the number of solvent molecules in the process of binding the ligand molecule to the target molecule, in accordance with the equation: Δribmd — -fcOomplex_Nsolventmol
"■" JN- lisolventmol
— V-tltarget_Lsolventmol ~*~ ■t'lιgand_Msolventιr»oI )
> wherein E
coπιpiex_:N_o.ve
nt
moi is the total energy of the target-ligand complex, bound to N molecules of the solvent, E
target_
Lsoι
ventπιoι is the total energy of the target molecule bound to L molecules of the solvent, E|
lgand_
Msoi
ventmoi is the total energy of the ligand molecule bound to M molecules of the solvent, E
SO|
V.
ntmoι is the energy of removal of one water molecule from the solvent, this energy calculated in the same quantum-mechanical approximation as the total energies of the intermolecular complex, target molecule and ligand molecule. In the case where balance of the solvent molecules is not observed, i.e., when N = L + M - K, in order to observe the energetic balance in the formula, in respect to which the enthalpy of binding the target molecule and the ligand olecule is calculated, the energy K*E
moi
S0|
veπt should also be taken into account. Calculation of the entropy of binding the ligand molecule to the target molecule in a solvent is carried out within the frame of the used quantum-mechanical approximation for those models of the intermolecular complex, the target molecule and the ligand molecule, using the same input files which were used in the calculations of the total energies of these molecular systems. This entropy of binding is determined by the vibrational frequencies of the molecules participating in the process under consideration. In order to calculate the vibrational frequencies in steps 64, 65 and 66 for the intermolecular complex, target molecule and ligand molecule, the force constants or second derivatives of the total energies of these molecular systems are calculated in respect to all the coordinates of the atoms in the points of corresponding minimums (global or local). A matrix (Hesse matrix) is formed from these values, and the vibrational frequencies are found by solving problems concerning the values of the Hesse matrix themselves. In accordance with the values of these frequencies, a calculation is made in step 67 (according to the equations, for example, in [10]) of the entropy component which is related to a loss of the vibrational degrees of freedom during the transition from the free states of the target molecule and the ligand molecule in the solvent to their intermolecular complex. The entropy of binding, due to losses of the vibrational degrees of freedom ΔSb,
n _v
ιbratιonai is calculated according to the equation:
), in the case of observance of the balance of the solvent molecules, and according to the formula: ^"bind _vιbratιonal ~ ϊ>complex_Nsolventmol
"■" Λ- ^solventmol
_ target solvent ol ' "lιgand_Msolventraol ),
in the case where the balance of the solvent molecules is not observed. A contribution to the change of entropy during the binding of the ligand molecule to the target molecule is provided by the fact that during binding these molecules each lose three forward and three rotating degrees of freedom, which they have in the free state. This part of the entropy is calculated in step 68 by the Sakur-Titrode equation, and only the coordinates and atomic weights of the atoms of the target molecule, ligand molecule and their intermolecular complex are used for this (see, for example [11]). Wherein the atomic weights of the atoms are unambiguously determined by their type (i.e., by their position in the D.I. Mendeleev Table), which is recorded in the initial files read from the information carrier in step 55 (Fig. 4). During the binding of the relatively small ligand molecules (usually these are molecules with an atomic weight that does not exceed 500 daltons) to target macromolecules such as proteins, the target molecule may be considered to be stationary, and the change of the entropy, calculated in step 6S, will be determined only by the characteristics of the ligand molecule. In step 69, using the values calculated in steps 67 and 68, the full change of entropy during the binding of the ligand molecule to the target molecule is calculated as the sum of changes of entropy obtained in steps 67 and 68. The process of binding the ligand molecule to the target molecule takes place in a solvent. Each of these molecules, being placed in the solvent, creates a cavity therein, pushing aside the solvent molecules. This process requires the expenditure of a certain free energy, which depends on the properties of the solvent. This free energy is called the cavitation energy, and it depends only on the properties of the solvent and the form of the molecule placed in the solvent. At present some methods are known for calculating the cavitation energy of solution (see, for example, the paper [12] and the references cited therein). Since during the binding of a ligand molecule to a target molecule and the formation of an intermolecular complex the form of the cavities in which the molecules under consideration are found changes, the cavitation free energy will also change. Its change ΔG
cav is calculated in step 70, and the result, being added in step 71 to the results of the calculations of the enthalpy of binding ΔH
bi
nd, obtained in step 63 and the entropy contribution - TΔSy
nd, calculated in step 69, provides the value of the total free energy of binding the ligand molecule to the target molecule in the solvent in accordance with the equation: ΔGbind
= ΔH iή - TΔSbind + ΔG
cav. In accordance with the instant invention, calculation of the enthalpy of binding the target molecule and the ligand molecule (step 63 in Fig. 4) may be carried out, using another method, in which the total energies of the target molecule, ligand molecule and their intermolecular complex are not calculated in the minimum point, but rather as the average values of the total energies obtained for a number of configurations of the intermolecular target-ligand complex.
Wherein the coordinates of the atoms of the target molecule, ligand molecule and their intermolecular complex, which correspond to each of such configurations, may be obtained within the frame of molecular-dynamic modeling. A block diagram is presented in Fig. 5 that illustrates a method of molecular-dynamic formation of the coordinates of atoms of the target molecule, ligand molecule and the intermolecular target-ligand complex for a set of configurations of the intermolecular complex, and also a method of calculating the enthalpy of binding the target molecule to the ligand molecule on the basis of average values of the total energies obtained for the set of those configurations. The initial data for the formation of input files for a set of configurations of an intermolecular complex, target molecule and ligand molecule are the coordinates, read in step 72 from the information carrier, of the atoms of the target molecule, ligand molecule and their intermolecular complex, which may be obtained from the PDB database or with the aid of different program builders of molecular structures. Furthermore, in the case where an explicit model of a solvent is used, the initial data are the coordinates of the atoms of the solvent molecule, which are introduced in step 73, in a manner similar to how that was carried out earlier in steps 27, 28 and 29. In the case where an implicit model of the solvent is used, which is introduced in step 74, similarly to the already described step 53, the coordinates of the atoms of the solvent molecules are not explicitly included in the modeling. The procedure of the molecular-dynamic formation of the input files with coordinates of the atoms of the target molecule, ligand molecule and intermolecular target-ligand complex for a number of configurations with an explicit model of the solvent includes the following steps. Using the coordinates of the atoms of the intermolecular complex of the target molecule and the ligand molecule, and also the solvent molecules (steps 26, 27 and 30 in Fig. 2) and, where this is necessary, introducing additional atoms to establish the necessary charge of the intermolecular target-ligand complex (steps 46, 49 in Fig. 3), molecular-dynamic modeling of the target-ligand complex is carried out in step 75 within the frame of the explicit model of the solvent. Wherein, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the target molecule and the ligand molecule and the presence of molecules of the solvent, the state of the molecular system is determined at which it will reach thermal balance. Then the paths of movement of all the atoms of the system are determined after they have reached thermal balance, when the average meanings of the values being observed (for example, the energy of interaction of the ligand molecule with the target molecule) stop changing in time (usually this is the length of the trajectory in time of about one nanosecond), and the coordinates of all the atoms of the interacting target molecule, ligand molecule and solvent molecules are recorded in step 76 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the
intermolecular target-ligand complex. Using the coordinates of the atoms of the target molecule, and also molecules of the solvent (steps 17, 28 and 31 in Fig. 2), and, where this is necessary, introducing additional atoms to establish the necessary charge of the target molecule (steps 47, 49 in Fig. 3), molecular- dynamic modeling of the target molecule is carried out in step 77 within the frame of an explicit model of the solvent. Vherein, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the target molecule and the presence of molecules of the solvent, the state of the molecular system is determined at which it will reach thermal balance. Then the paths of movement of all the atoms of the system are determined after they have reached thermal balance and the coordinates of all the atoms of the interacting target molecule and solvent molecules are recorded in step 78 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the target molecule. Using the coordinates of the atoms of the ligand molecule, and also molecules of the solvent (steps 24, 29 and 32 in Fig. 2), and, where this is necessary, introducing additional atoms to establish the necessary charge of the ligand molecule (steps 48, 49 in Fig. 3), olecular- dynamic modeling is carried out in step 79 within the frame of an explicit model of the solvent. Wherein, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the ligand molecule and the presence of molecules of the solvent, the state of the molecular system is determined in step 79 at which it will reach thermal balance. Then the paths of movement of all the atoms of the system are determined after they have reached thermal balance and the coordinates of all the atoms of the interacting ligand molecule and solvent molecules are recorded in step 80 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the ligand molecule. The procedure of the molecular-dynamic formation of the input files with coordinates of the atoms of the target molecule, ligand molecule and intermolecular target-ligand complex for a number of configurations with an implicit model of the solvent includes the following steps. Using the coordinates of the atoms of the intermolecular complex of the target molecule and the ligand molecule (step 26 in Fig. 2) and, where this is necessary, introducing additional atoms to establish the necessary charge of the intermolecular target-ligand complex (steps 46, 49 in Fig. 3), molecular-dynamic modeling is carried out within the frame of the implicit model of the solvent. The solvent in this case is taken into account with the aid of an implicit model as described above (step 53 in Fig. 3) within the frame of solution of the Poisson, Poisson- Boltzmann equations, the equation for continuous conductor or by another method, talcing into account the effect of the solvent without using the coordinates of the solvent molecules. Then, in
a manner similar to how this was done within the frame of the explicit model of the solvent, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the target molecule and the ligand molecule and the presence of molecules of the solvent, the state of the molecular system at which it will reach thermal balance is determined in step 75. Then the paths of movement of all the atoms of the system are determined after they have reached thermal balance and the coordinates of all the atoms of the interacting target molecule, ligand molecule and solvent molecules are recorded in step 76 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the intermolecular target-ligand complex. Using the coordinates of the atoms of the target molecule and, where this is necessary, introducing additional atoms to establish the necessary charge of the target molecule (steps 47, 49 in Fig. 3), molecular-dynamic modeling is carried out within the frame of the implicit model of the solvent. The solvent in this case is taken into account with the aid of an implicit model (step 53 in Fig. 3) within the frame of solution of the Poisson, Poisson-Boltzmann equations, the equation for continuous conductor or by another method, taking into account the effect of the solvent without using the coordinates of the solvent molecules. Then, in a manner similar to how this was done within the frame of the explicit model of the solvent, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the target molecule and the presence of molecules of the solvent, the state of the molecular system at which it will reach thermal balance is determined in step 77. Then the paths of movement of all the atoms of the system are determined after they have reached thermal balance and the coordinates of all the atoms of the interacting target molecule and solvent molecules are recorded in step 78 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the target molecule. Using the coordinates of the atoms of the ligand molecule and, where this is necessary, introducing additional atoms to establish the necessary charge of the ligand molecule (steps 47, 49 in Fig. 3), molecular-dynamic modeling is carried out within the frame of the implicit model of the solvent. The solvent in this case is taken into account with the aid of an implicit model (step 53 in Fig. 3) within the frame of solution of the Poisson, Poisson-Boltzmann equations, the equation for continuous conductor or by another method, taking into account the effect of the solvent without using the coordinates of the solvent molecules. Then, in a manner similar to how this was done within the frame of the explicit model of the solvent, by solving the Newton or Lagrange equations, taking into account the thermal movement of all the mutually interacting atoms of the ligand molecule and the presence of molecules of the solvent, the state of the molecular system at which it will reach thermal balance is determined in step 79. Then the paths of movement of all the atoms of the system are determined after they have reached thermal
balance and the coordinates of all the atoms of the interacting ligand molecule and solvent molecules are recorded in step 80 on the information carrier after certain periods of time a certain number of times for a predetermined set of configurations of the ligand molecule. Further, for each set of coordinates of the intermolecular target-ligand complex in the explicit or implicit model of the solvent obtained in step 76, all the atoms are divided in step 81 into groups in a manner similar to the way it was earlier done (step 35 in Fig. 35). For each set of coordinates of the target molecule in the solvent, obtained in step 78, all the atoms are divided into groups in accordance with method 82, similar to method 36. For each set of coordinates of the ligand molecule in the solvent obtained according to method 80, all the atoms are divided in step 83 into groups in a manner similar to the way it was earlier done (step 37 in Fig. 35). For each configuration of the intermolecular target-ligand complex, i.e., for each set of its coordinates which is obtained in molecular-dynamic modeling, an input file for quantum- mechanical modeling of the intermolecular target-ligand complex in the solvent is formed in step 84, as described above (step 4 in Fig. 1 and step 56 in Fig. 4). Similarly, for each configuration of the target molecule, i.e., for each set of its coordinates which is obtained in molecular-dynamic modeling, an input file for quantum- mechanical modeling of the target molecule in the solvent is formed in step 85, as described above (step 4 in Fig. 1 and step 56 in Fig. 4). For each configuration of the ligand molecule, i.e., for each set of its coordinates which is obtained in molecular-dynamic modeling, an input file for quantum- mechanical modeling of the ligand molecule in the solvent is formed in step 86, as described above (step 4 in Fig. 1 and step 56 in Fig. 4). For each input file formed in step 84, i.e., for each configuration of an intermolecular target-ligand complex, with the aid of the method of modeling selected in step 87, calculation of the total energy of the intermolecular target-ligand complex in the solvent is carried out in step
88 without optimization of the spatial structure of the complex. Then in step 89 the average value of the total energy of the intermolecular target-ligand complex in the solvent is calculated in respect to the whole selected set of configurations of that molecular system. For each input file formed in step 85, i.e., for each configuration of a target molecule, with the aid of the method of modeling selected in step 87, calculation of the total energy of the target molecule in the solvent is carried out in step 90 without optimization of its spatial structure. Then in step 91 the average value of the total energy of the target molecule in the solvent is calculated in respect to the whole selected set of configurations of that molecular system. For each input file formed in step 86, i.e., for each configuration of a ligand molecule,
with the aid of the method of modeling selected in step 87, calculation of the total energy of the ligand molecule in the solvent is carried out in step 92 without optimization of its spatial structure. Then in step 93 the average value of the total energy of the ligand molecule in the solvent is calculated in respect to the whole selected set of configurations of that molecular system. Using the obtained average values of the total energies of the intermolecular target- ligand complex in the solvent, the target molecule in the solvent and the ligand molecule in the solvent, obtained in steps 89, 91 and 93, respectively, in step 94 the enthalpy of binding the ligand molecule to the target molecule in the solvent is calculated as the difference between the average value of the total energy of the intermolecular target-ligand complex in the solvent, on the one hand, and the sum of the average values of the total energies of the target molecule in the solvent and the ligand molecule in the solvent, on the other hand. Further, the entropy of binding the ligand molecule to the target molecule in the solvent and the free energy of binding the ligand molecule to the target molecule in the solvent are calculated (steps 69 and 71 in Fig. 4) . The method presented above for modeling the binding of ligand molecules to target molecules by methods of quantum mechanics, taking the effect of the solvent into account, was carried out by the authors of the instant invention for a number of intermolecular complexes of proteins and ligands, for which the energetic characteristics of binding were experimentally known, in particular, the enthalpy of binding. In accordance with the claimed method, the calculation of the enthalpy of binding ligand proteins was carried out for these test complexes. By comparing the obtained results with corresponding experimental values, an initial assessment of the accuracy of the claimed method for modeling was provided. In order to carry out test modeling, eight complexes of proteins with ligands were selected, in respect to which not only the experimental values of the enthalpy of binding were known, but also the structures obtained with the aid of X-ray structural analysis and entered into the PDB (Protein Data Bank) database [7]. A list of studied complexes is presented in Table 1 . In this table, the names of these complexes in accordance with the PDB database, the Latin names of corresponding proteins and ligands, the experimental values of the enthalpy of binding for each complex, and also citations to papers from which these data were obtained, are presented. The structural formulas and structures of ligands of corresponding complexes optimized within the frame of quantum mechanics are provided in Fig. 6.
Table 1 Complexes of ligand proteins with experimentally known values of the enthalpy of binding, which were selected as test systems for modeling binding in the ligand-protein system by methods of quantum mechanics
As is evident from the presented Table 1 and from Fig. 6, the selected complexes comprise ligands that are different organic compounds in respect to chemical makeup and properties. Among them are carbohydrates (systems 1, 2, 6 and 7 in the numbering of Table 1), heterocyclic compounds (systems 3 and 4 in the numbering of Table 1), an aromatic compound (system 5 in the numbering of Table 1) and an oligopeptide (system 8 in the numbering of Table 1). Some of these ligands are negatively charged, some are neutral, the number of atoms therein vary from 20 atoms to 94. The experimental values of the enthalpies of binding change in the selected row of complexes from -2.7 to -28.6 kcal/mol. Two selected complexes (systems 3 and 4 according to the numbering of Table 1) comprise identical ligands and different proteins, two
other complexes (systems 6 and 7 according to the numbering of Table 1) comprise different ligands and identical proteins. Thus, the selected series of intermolecular ligand-protein complexes was a sufficiently wide series of systems, which is of interest during the testing of the claimed method of modeling. The initial structures of the complexes, in particular, the Cartesian coordinates of their atoms were taken from the PDB database. For all practical purposes, all the initial structures of the proteins and ligands contained only coordinates of heavy (C, N, O, S) atoms. Hydrogen atoms were added to the structures of the proteins of all the complexes with the aid of a REDUCE program [19] in the first step of building models of the complexes. Hydrogen atoms were added to the structures of the ligands of all the complexes with the aid of the Materials Studio Visualizer packet, provided by the Accelrys company [20]. In this case, taking the effect of the solvent into account was carried out within the frame of the explicit model of water, but in the instant modeling the addition of water molecules in molecular mechanics or in molecular dynamics was not carried out. The positions of water molecules, bound to the complex and its components, were taken from experimental data (the PDB data of X-ray-structural analysis) on the structured water molecules in the complex. Sorting the atoms into two groups in accordance with the method proposed in the claimed invention was carried out for each complex. This procedure was carried out with the aid of the LigEnv program [21], developed in the Algodign company. Ligand atoms and protein atoms were included in the first group for each complex, these atoms composing whole amino acid residues, at least one of which being at a distance, not exceeding a predetermined value, from each ligand atom. Furthermore, the first group included atoms of water molecules, which atoms were present in the input data obtained from the PDB database, each of which being at a distance from the ligand atoms that did not exceed a predetermined value. All the remaining atoms of the protein and water molecules were included in the second group. During the building of the protein model, protein atoms, relating to the first group of atoms for the intermolecular complex, were included in the first group of atoms. In order to model the enthalpy of binding in all the selected complexes, an approach was selected in which the atoms of the protein, ligand and their complex together with water molecules, relating to the first groups, were modeled by methods of quantum mechanics, while the atoms of the second groups were withdrawn from modeling. With such a method of modeling, as already indicated, in the first place it was necessary to select sufficiently large protein models, which adequately describe the place of binding the ligand, and, consequently, in order to accelerate the count, a semiempirical quantum- mechanical method of calculating PM3 [22] was selected. This method, together with a sufficiently high speed of counting provides satisfactory accuracy during calculation of the total energy of molecular systems that have hydrogen bonds. In the second place, with such an
approach, it was necessary to introduce a number of additional atoms, which restored disrupted covalent bonds of protein atoms during their division into groups. The disrupted covalent bonds of the N-ends of the main chains of proteins were joined by -C(0)-(CH3) groups, while the disrupted covalent bonds of the C-ends of the main chains of proteins were joined by — (H)- (CH3) groups. The additionally added atoms whereby were introduced together with the atoms of the first group into the modeling by methods of quantum mechanics. The coordinates of atoms of ligands, water molecules and side chains of proteins were selected as the parameters to be optimized. In order to maintain the three-dimensional structure of the protein, which structure was taken from the experiment, the coordinates of heavy atoms (C, N, O) of the main chains of the proteins and the coordinates of the added atoms were presumed to be non-optimizable parameters. The charge states of different molecular groups of proteins and ligands in water with known pH values were determined in the course of forming the input data for modeling. In accordance with this, additional protons in corresponding acidic or basic amino groups of the protein were added or withdrawn. By summing the charge states of the acidic and basic amino groups, the net charge of the models of the proteins, ligands and their intermolecular complexes was determined. Fig. 7 shows an optimized structure of a model of one of the studied complexes built in such a manner. The structures of all the remaining models, in the general case, look similar and are left out of consideration. In the course of modeling, it was shown that special attention should be given to modeling water as a solvent taken into account in an explicit manner. Two methods were proposed for taking into account the disruption and formation of hydrogen bonds during the transition from the complex to its components and the participation of water molecules in that process. Accordingly, two methods of adding water molecules during the transition from a complex to its components were given consideration. Conditionally, these methods were called approaches WI and W2. In accordance with approach WI, water molecules, the position of which was taken from an experiment from the PDB database, bound in a complex to a protein, were counted during the modeling of individual components of the complex together with the protein. Water molecules, the position of which was taken from the experiment from the PDB database, bound in the complex to a ligand, were counted during the modeling of individual components of the complex together with the ligand. Water molecules, bound in the complex simultaneously to a protein and to a ligand and forming with them bridges with hydrogen bonds, were taken into account during calculation of the protein and ligand separately two times. In accordance with approach W2, water molecules, the position of which was taken from
an experiment, were also taken into account by a method similar to the method WI. Furthermore, in that case additional molecules of water were added to the protein and to the ligand in those places where there were hydrogen bonds in the complex between the protein and the ligand, so that disrupted hydrogen bonds are replaced during a transition from the intermolecular complex to its components. A schematic presentation of these two methods WI and W2 for taking account of water molecules in an explicit manner is shown in Fig. 8 and Fig. 9, respectively. In accordance with the presented schemes, calculation of the enthalpy of binding a protein and ligand in a complex in accordance with method WI was carried out in accordance with the equation: ΔHbind (WI) = ΔHLp - (ΔHp[nW(P), AW( P)J + ΔHL[ΠIW( ), *W(LP)]) + k ' ΔH • Here ΔHLp is the calculated energy of the ligand-protein (LP) complex taken together with all the water molecules interacting with the complex. ΔHP[n (p), AW(LP)] is the calculated energy of the protein (P) taken together with all the (n + k) water molecules interacting therewith (see Fig. 8), ΔHL[mW(L), A ( P)]) is the calculated energy of the ligand (L), taken together with all the (m + k) water molecules interacting therewith (see Fig. 8). k " ΔHW is the calculated energy of k water molecules W(LP) (see Fig. 8), which should be added to the sum to maintain energetic balance. ΔHw is the energy of a water molecule, which is calculated, taking into account its withdrawal from a drop of water. In this case the calculated energy ΔHw = -53.43 kcal/mol. Calculation of the enthalpy of binding the protein and ligand in a complex in accordance with the method W2 was carried out in accordance with the equation: ΔHbind (W2) = ΔH P - (ΔHp[nW(P), *W(LP), qW(PHB)] + ΔHL[mW( ), AW(LP), /W(LHB)]) + (k+q+t) ' ΔHw • Here ΔHL is the calculated energy of the ligand-protein (LP) complex taken together with all the water molecules interacting with the complex. ΔHP[nw(p), AW( P), W(PHB)] is the calculated energy of the protein (P) taken together with all the (n+k+q) water molecules interacting therewith (see Fig. 9), ΔH [ΠIW(L), A ( P), .W(LHB)]) is the calculated energy of the ligand (L), taken together with all the (m+k+t) water molecules interacting therewith (see Fig. 9). (k+q+t) ' ΔHw is the calculated energy of (k+-q+t) water molecules W(LP), W(PHB) and W(LHB) (see Fig. 9), which should be added to the sum to maintain energetic balance. The enthalpies of binding proteins and ligands for complexes taken as test systems (see
Table 1, systems 1-8), obtained in approximations WI and W2 in accordance with the claimed invention, are presented in Table 2 together with corresponding experimental data and deviations, calculated from experimental values. As is evident from that table, in both methods of taking account of the water in an explicit manner, a sufficiently sensible agreement was obtained between the calculated and
experimentally obtained data in respect to the enthalpy of binding in ligand-protein complexes. The root-mean-square (RMS) deviation of the enthalpy of binding calculated for all the complexes in the approach WI was 7.0 kcal/mol, a similar root-mean-square deviation obtained for the method W2 was 1.85 kcal/mol. The latter accuracy of calculation of the enthalpy of binding is sufficiently high and comparable with the accuracy of existing experimental method of measuring the enthalpy of binding in ligand-protein complexes. Table 2 Calculated enthalpy of binding proteins and ligands for complexes 1-8 in accordance with Table 1
122, 534-535. 3. Lyne, P.D.; Hodoscek, J.; Karplus, M., A hybrid QM/MM potential employing Hartree-Fock of density functional methods in the quantum region, J. Phys. Chem. A., 1999, 103, 3462-3471. 4. Lewis, J.P.; Carter, C.W., Jr.; Hermans, J.; Pan, W.; Lee, T.-S.; Yand, W., Active species for the ground-state complex of cytidine deaminase: a linear-scaling quantum mechanical investigation, J. Am. Chem. Soc, 1998, 120, 5407-5410. 5. York, D.M.; Lee, T.-S.; Yang, W., Quantum mechanical treatment of biological macromolecules in solution using linear-scaling electronic structure methods, Phys. Rev. Lett., 1998, 80, 5011-5014. 6. Lee, T.-S.; York, D.M.; Yang, W., Linear-scaling semiempirical quantum calculations for macromolecules, J. Chem. Phys., 1996, 105, 2744-2750. 7. Berman, H.M.; Westbroo , J.; Feng, Z.; Gilliland, G.; Bhat, T.N.; Weissig, H.; Shindyalov, I.N.; Bourne, P.E., Nucleic Acids Research, 28 (2000), 235-242, http://www.rcsb.org/pdb. 8. Murphy, R.B.; Philipp, D.M.; Friesner, R.A., J. Comp. Chem., 2000, 21, 1442-1457. 9. Stewart, J.J.P., J. Mol. Struct. (Theochem), 1997, 401, 195-205. 10. Tidor, B.; Karplus, M., The contribution of vibrational entropy to molecular association, J.Mol.Bio., 238 (1994) 405-414. 11. Sonja Schwartzl, et al., J. Comput. Chem.., 23 (2002), 1143-1149. 12. Levy, R.M.; Gallicchio, E., Ann. Rev. Phys. Chem., 49 (1998), 531. 13. Fukada, H.; Sturtevant, J.M.; Quiocho, F.A., J. Biol. Chem., 1983, 258, 13193- 13198. 14. Schwarz, F.P.; Ahmed, H.; Bianchet, M.A.; Amzel, L.M.; Vasta, G.R., Biochem, 1998, 37, 5867-5877. 15. Hyre, D.E.; Trong, I.L.; Freitag, S.; Stenkamp, R.E.; Stayton, P.S., Prot. Sci, 2000, 9, 878-885. 16. Holdgate, G.A.; Tunnicliffe, A.; Ward, W.H.; Weston, S.A.; Rosenbrock, G.; Barth, P.T.; Taylor, I.W.; Paupti,t R.A.; Timms, D., Biochem. 1997, 36, 9663-9673. 17. Klotz, I.M., Ligand-Receptor Energetics/ a guide for the perplexed, John Wiley &
Sons, New York, Chichester, Weinheim, Brisbane, Singapore, Toronto, 1997, p. 77.
18. Henriques, D.A.; Ladbury, J.E.; Jackson, R.M., Prot. Sci. 2000, 9, 1975-1985. 19. Word, J.M.; J. Mol. BioL, 1999, 285, 1733-1745. 20. Materials Studio Visualizer, Materials Studio Release 2.2, San Diego: Accelrys Inc., 2002, http://www.accelrys.com/mstudio/visualizer.html 21. Vlasov, P., LigEnv: program for construction of a model represented an active site of protein-ligand complex from PDB molecular structure, Algodign, LLC, Moscow, 2002. 22. Stewart, J.J.P, J. Comput. Chem.., 1989, 10, 2, 209-264.