# Patent

## Title
Method for searching heterogeneous compound databases using topomeric shape descriptors and pharmacophoric features

## Abstract

Heterogeneous compound databases can be searched for compounds which are likely to have the same biological activity as a known (query) molecule. Query molecules and the molecules in the database are split into fragments according to common fragmentation rules. Fragments are aligned in a uniform conformation according to a topomeric alignment process and interaction energy fields, typically steric fields, between a probe and the fragment atoms are generated to capture the fragment shapes. Comparison of the fields for the query fragments with the fields for the database compound fragments yields a measure of shape similarity. Searches for similarly shaped substructures and cores can also be readily accomplished. Pharmacophoric style features can be defined for the topomerically aligned fragments but with user specified weighting of the importance of each. Differences in features are defined with the same dimensionality as shape so that both shape and features can be used to search.

## Background
A portion of the disclosure of this patent document contains material which is subject to copyright protection. The copyright owner has no objection to the facsimile reproduction by anyone of the patent document or the patent disclosure, as it appears in the Patent and Trademark Office patent file or records, but otherwise reserves all copyright rights whatsoever.
BACKGROUND OF THE INVENTION
1. Field of the Invention
This invention relates generally to the field of pharmaceutical research and to the three dimensional searching of structures of chemical compounds to identify compounds which may share a biological activity with a known compound. In particular the invention concerns a method for searching databases of commercially available compounds which may or may not share any common synthetic linage.
2. Description of Related Art
The advent of high throughput screening of chemical compounds for biological activity has dramatically changed the paradigm of pharmaceutical research in recent years. Coupled with combinatorial synthesis, it is now possible to test millions of compounds on an efficient basis. However, the cost per hit of such searching remains extremely high given the enormous number of compounds which can be tested and the typically low “hit” rates which are achieved. As a result, greater emphasis has been placed on the testing of compound libraries which are believed to contain a higher percentage of potentially relevant molecules. The skills of computational chemists have been employed to design such compound libraries for testing.
Two type of libraries were considered possible: first, a library which explored the diversity of structures in chemical space across the range of compounds which could be synthesized without oversampling the same area of diversity space (redundant testing); and second, a library in which the compounds would be likely to have the same biological activity as a known molecule or drug. The major problem confronting computational chemists in the selection of compounds for such libraries was how to characterize the compounds in a manner which would permit the desired selections. Bioscientists have long known that the three dimensional shape of a compound which acts as a ligand to a larger biomolecule must be complimentary to the shape of the binding site of the larger biomolecule. In studying the relationships between the chemical structure of a molecule and its biological activity (structure activity relationships [SAR]) many techniques to characterize the three dimensional shape of molecules were devised. One of the most successful of the techniques for generating a quantitative structure activity relationship (QSAR) characterized the shape of molecules by defining an interaction energy field between a probe molecule and each part of the studied molecule in a three dimensional grid surrounding the molecule. The shape data thus generated for a series of molecules could be correlated with the biological activity of the molecules to produce the QSAR. This technique by Cramer and Wold (Comparative Molecular Field Analysis [CoMFA]) is described in detail in U.S. Pat. No. 5,025,388 and U.S. Pat. No. 5,307,287.
Use of the CoMFA approach required detailed considerations of two major factors: 1) the proper alignment of the test molecules; and 2) the conformation or conformations of the molecules which had to be taken into account. In addition, the technique worked only with molecules sharing the same biological activity. However, the technique clearly demonstrated the power of utilizing three dimensional shape descriptors in molecular analysis.
Over time many three dimensional shape descriptors and methods of library selection were attempted by computational chemists. U.S. Pat. No. 5,703,792 to Chapman describes one such approach. Two major problems confronted the field and cast doubt on the generality or accuracy of all the methods which had been devised. The first problem was that no one could show that the molecular structural descriptors which had been used were generally valid; that is, that the descriptors described molecules in a manner which correlated with biological activity across a range of biological systems. Any descriptor which would be used to select compounds for libraries would have to be valid irrespective of the biological activity which might be tested against the library. The second problem was that there was likewise no way to demonstrate that the methods of handling multiple conformations in the prior art methods were either accurate or applicable across all types of molecules.
The solution to these problems by Cramer, Patterson, Clark, and Ferguson are taught in U.S. Pat. No. 6,185,506. The validity of a molecular structural descriptor can be demonstrated across multiple biological activities by employing the Patterson plot methodology described in the patent. Both two and three dimensional descriptors can be evaluated by the methodology, and, in principal, there is no limitation on the dimensionality of the descriptors which can be evaluated. Using the validation technique, valid descriptors were identified which could be used with assurance to design libraries having desired properties. By this method the two dimensional prior art fingerprint Tanimoto descriptor was shown to be valid as well as a new three dimensional descriptor described below. The validation methodology also identified a neighborhood distance characteristic of the descriptors which could be used in the design of the libraries. In addition, the neighborhood distance led directly to methods for searching the libraries, and, once a molecule had shown activity in a screen, for expanding the search for other molecules having the same activity.
Further, a solution to the problem of identifying a generally appropriate molecular conformation or conformations to take into account was taught. An alignment rule for molecular parts (topomeric alignment) is demonstrated which generates a uniform orientation. The shape of the molecular part is characterized, as in CoMFA, by a field of interaction energies calculated between a probe and the atoms in the aligned molecular part at each point in a three dimensional grid surrounding the molecular part. The steric interaction energies are principally used although, in the appropriate circumstances, electrostatic interaction energies may be added. Although the alignment may be arbitrary and unlikely for any particular molecule, the field shape descriptor of the topomeric alignments was shown to be a valid molecular structural descriptor by means of the Patterson plot method.
Using descriptors having an associated neighborhood distance, molecules could be identified which shared shape characteristics in a way which was meaningfully related to their biological activity. The problems of efficient library design and selection of combinatorially accessible molecules could be further addressed. In U.S. patent application Ser. No. 08/903,217, presently allowed, the construction and searching of a virtual library is described. The virtual library contains validated molecular structural descriptions of each component part which could be used in a specified combinatorial synthesis. All possible product molecules which could be combinatorially derived from the component parts can be searched, without the necessity of generating the product structures during the search, for product molecules having desired properties by searching through only a combination of the descriptors of the component parts of the product molecules. In the preferred embodiment the Tanimoto and the three dimensional topomeric CoMFA descriptors are employed.
Due to the combinatorial nature of the number of product molecules whose characteristics can be determined, a relatively small number of structural variations (tens of thousands), cores, and synthetic schemes employing only two attachment points can yield a searchable library of billions of possible molecules according to the method of the patent. Indeed, the number of searchable molecules outnumbers the number of molecules ever reported by several orders of magnitude. By the techniques disclosed in the patent, this virtual library can be searched very fast to construct diverse libraries of molecules likely to share the same biological activity or to find molecules which share the same biological activity as a combinatorially derived query molecule. Further, query molecules which derive from unknown synthetic routes can be fragmented and the molecular descriptor characterization of the fragments used to search for similarly shaped fragments and potential molecules with likely similar biological activity defined in the virtual library. In practice the topomeric field molecular structural descriptor has proven to be very valuable in searching the virtual library. The powerful and fast searching capabilities of the virtual library method have yielded significant advances.
However, the molecules in the virtual library which can be searched by definition derive from a combinatorial assembly of a relatively few number of constituent parts and can be said to be homogeneous in that sense. By virtue of the exceedingly large size of the virtual library, 5 molecules may be identified which are not readily available. Also, although the possible product molecules which can be searched are the result of known combinatorial synthetic schemes, the actual synthesis may not be easily achieved. In the day to day world of pharmaceutical research, large assemblages of available molecules can be commercially obtained. These assemblages are not the result of any particular combinatorial synthesis but rather represent the assembly of a wide range of molecules from many different sources and syntheses, some known, some unknown. Therefore, these assemblages of molecules can be characterized as heterogeneous.
It would be useful if heterogeneous assemblages of available molecules could be searched for molecules which are likely to have a biological activity similar to a known compound before synthesis of new compounds is undertaken with the concomitant additional time and expense.
BRIEF SUMMARY OF THE INVENTION
Databases which contain the structures of a heterogenous assembly of available molecules can be searched for molecules having a biological activity similar to a known compound. Each molecule specified by the database is split into several fragments according to defined rules and the shape of those fragments is compared to the shape of the fragments generated from a query molecule using the topomeric field molecular structural descriptor. The molecules having the closest matching shapes to the query molecule are selected for further testing.



BRIEF DESCRIPTION OF THE FIGURES
FIG. 1 shows a number of possible ways to fragment a molecule into two pieces in accordance with the fragmentation rule.
FIG. 2 shows a number of possible ways to fragment a molecule into three pieces in accordance with the fragmentation rule.



DETAILED DESCRIPTION OF THE INVENTION
Computational Environment:
Generally, all calculations and analyses to perform the method of the disclosed invention are implemented in a modern computational chemistry environment using software designed to handle molecular structures and associated properties and operations. For purposes of this Application, such an environment is specifically referenced. In particular, the computational environment and capabilities of the SYBYL, UNITY, and CONCORD software programs developed and/or marketed by Tripos, Inc. (St. Louis, Mo.) are specifically utilized. The software code to implement the method of the disclosed invention is set out in the Appendices to this Application. Software with similar functionalities to SYBYL, UNITY, and CONCORD are available from other sources, both commercial and non-commercial, well known to those in the art. A general purpose programmable digital computer with ample amounts of memory and hard disk storage is required for the implementation of this invention. In performing the methods of this invention, representations of thousands of molecules and molecular structures as well as other data may need to be stored simultaneously in the random access memory of the computer or in rapidly available permanent storage. The inventors use Silicon Graphics, Inc. (SGI) “R12000” computers having 350-400 MHz processors and between 256 Mb and 512 Mb of memory with 8-10 Gb hard drive storage disks. In addition SGI “Origin” or “O2” or “O2100” computers can be used. Access to several gigabytes of storage and faster Silicon Graphics, Inc. processors is useful.
Incorporation of Patent Disclosures:
The disclosures of U.S. Pat. No. 6,185,506 and of U.S. patent application Ser. No. 08/903,217 are expressly and completely incorporated into this application as if fully set forth herein.
Topomeric Alignment:
As taught in the incorporated U.S. Patent and patent application, molecular fragments may be aligned following topologically-based rules to generate a single, consistent, unambiguous, aligned topomeric conformation. The procedure also takes full account of chiral atoms. All fragments which are to be compared in a search must be aligned with the same topomeric rules. In the present method such a topomeric alignment is used, the details of which are fully set out in the attached software code.
Calculation of Fields:
The basic CoMFA methodology provides for the calculation of both steric and electrostatic fields. It has been found up to the present point in time that using only the steric fields yields a better molecular structural descriptor than a combination of steric and electrostatic fields. There appear to be three factors responsible for this observation. First is the fact that steric interactions—classical bioisosterism—are certainly the best defined and probably the most important of the selective non-covalent interactions responsible for biological activity. Second, adding the electrostatic interaction energies may not add much more information since the differences in electrostatic fields are not independent of the differences in steric fields. Third, the addition of the electrostatic fields will halve the contribution of the steric field to the differences between one shape and another. This will dilute out the steric contribution and also dilute the neighborhood property. Clearly, reducing the importance of a primary descriptor is not a way to increase accuracy. However, it is certainly possible that in a given special situation the electrostatic contribution might contribute significantly to the overall “shape”. Under these unique circumstances, it would be appropriate to also use the electrostatic interaction energies or other molecular characterizers, and such are considered within the scope of this disclosure. In particular, as will be discussed below, it has been found that the additional information typically associated with pharmacophore mapping can be utilized to further characterize the similarity between topomerically aligned molecular fragments.
The steric fields of the topomerically aligned molecular fragments are generated almost exactly as in a standard CoMFA analysis using an sp3 carbon atom as the probe. In standard CoMFA, both the grid spacing and the size of the lattice space for which data points are calculated will depend on the size of the molecule and the resolution desired. Typically, a 2 Å grid spacing in employed both in CoMFA and in the heterogenous database searching method of the present disclosure. However the grid dimensions are varied in the present invention. For query molecules, the size of the grid is adjusted to encompass the smallest region that all of the query fragments will fit into. This significantly reduces the number of calculations that are necessary without reducing the ability of the descriptor to fully characterize the structures. This modification will be discussed in more detail below. The steric fields are set at a cutoff value (maximum value) as in standard CoMFA for lattice points whose total steric interaction with any side-chain atom(s) is greater than the cutoff value.
One difference from the usual CoMFA procedure is that atoms which are separated by one or more rotatable bonds are set to make reduced contributions to the overall steric field. An attenuation factor, preferably about 0.85, is applied to the steric field contributions which result from these atoms. For atoms at the end of a long molecule, the attenuation factor produces very small field contributions (ie: [0.85]N) where N is the number of rotatable bonds. This attenuation factor is applied in recognition of the fact that the rotation of the atoms provides for a flexibility of the molecule which permits the parts of the molecule furthest away from the point of attachment to assume whatever orientation may be imposed by the unknown receptor. If such atoms were weighted equally, the contributions to the fields of the significant steric differences due to the more anchored atoms (whose disposition in the volume defined by the receptor site is most critical) would be overshadowed by the effects of these flexible atoms.
Topomer Similarity:
The notion of topomer similarity between a pair of molecules is defined as the “distance” represented by the difference between the molecular fields which serve to characterize the molecules' shapes. As an example, assume two molecules A and B which have each been placed in their topomeric alignment and the steric field values calculated for each point in the surrounding three dimensional grids. Let each grid point be denoted by its corresponding Cartesian X, Y, Z coordinate so that for each molecule the grid points are defined as X0, Y0, Z0 . . . XN, YN, ZN.
For each molecule A and B the field values, VA and VB, at each point in the grid are denoted as:


    
    
        VAX0, VAY0, VAZ0 . . . VAXN, VAYN, VAZN and VBX0, VBY0, VBZ0 . . . VBXN, VBYN, VBZN.

The root sum square of distances between the fields is then defined as:

    
    






  
    
      
        
          
            
              (
              
                
                  V
                  X0
                  A
                
                :
                
                  V
                  X0
                  B
                
              
              )
            
            2
          
          +
          
            
              (
              
                
                  V
                  Y0
                  A
                
                :
                
                  V
                  Y0
                  B
                
              
              )
            
            2
          
          +
          
            
              (
              
                
                  V
                  Z0
                  A
                
                :
                
                  V
                  Z0
                  B
                
              
              )
            
            2
          
          +
          ...
          +
        
      
    
    
      
        
          
            
              (
              
                
                  V
                  XN
                  A
                
                :
                
                  V
                  XN
                  B
                
              
              )
            
            2
          
          +
          
            
              (
              
                
                  V
                  YN
                  A
                
                :
                
                  V
                  YN
                  B
                
              
              )
            
            2
          
          +
          
            
              (
              
                
                  V
                  ZN
                  A
                
                :
                
                  V
                  ZN
                  B
                
              
              )
            
            2
          
        
      
    
  




This distance is conveniently denoted as:

√{square root over ((A:B)2)}

For identical molecular structures, the distance equals 0. Therefore, the closer the value of the distance is to zero, the closer in shape two molecules will be. When searching among many possible structures, the minimum calculated value of the distance is sought.

Fragmentation:

The following critical question which frequently occurs in chemical research, and especially in biological research, can now be addressed. The problem, as it is usually presented, takes the form: given an arbitrary query molecule (generally one previously found to exhibit a desired activity), find biologically similar molecules, that is molecules of similar 3D shape and activity. Generally, such a query molecule will not have resulted from a combinatorial synthesis, and, in fact, no knowledge of a possible synthetic route to the molecule may be available. In searching the virtual library of application Ser. No. 08/903,217, the topomeric 3D shape data within the virtual libraries actually describe fragments (structural variations) of molecules. To find similarly shaped molecules within the virtual library, the query molecule must be fragmented and the shapes of its fragments compared with the shapes of corresponding fragments (structural variations) in the virtual library. The difficulty is that a query molecule can be fragmented in so very many ways. The solution adopted for virtual library searching was a way to emphasize those fragmentations that are most likely to conform to efficient synthetic routes from available starting materials, without requiring the searcher of the virtual library to have any knowledge of what synthetic routes it includes.
The solution employed a “fragmentation table”, where each row constitutes a rule of the following sort: “for each occurrence of this particular structural feature combination (structural variation) in the query molecule, decompose the query molecule in a particular way specified in terms of this structural feature, and search only those combinatorial libraries that utilize specified reactions (sequences) and/or building blocks, mapping specified query fragments onto specified classes of building blocks”. Each such query decomposition found generates a search of the virtual library, returning all those products whose sum of squares of differences in shape between corresponding product and query fragments is less than a user specified neighborhood distance threshold. Passing the query molecule (by means of a suitable computer program) against all the rows of this table generates all searches.
The situation is much more complicated when a search of a database of heterogeneous compounds is desired. Not only is it necessary to fragment the query molecule, but each molecule in the database has to be likewise fragmented and comparisons made between the query fragments and the fragments arising from each molecule. Typically, anywhere from 2 to 50 different fragments might be generated by fragmenting each molecule in the database. To compare 6 fragments from a query molecule to an average of 20 fragments from each of 50,000 molecules in a heterogeneous database would require 6×20×50,000=6,000,000 field comparisons. [Actually, as will be described below, because fragment pairs or triplets are involved, cross comparisons increase this number.] This is at least an order of magnitude greater than the typical 6 fragment query comparison to even 50,000 structural variations in the virtual library. In principal, a virtual library of every fragment occurring in all of the molecules in all examined heterogenous databases could be assembled, but the size of such a virtual library and the complexities of searching are not trivial.
The method adopted for the present invention does not precalculate and store the metric characteristics of each fragment of each heterogenous database molecule. Rather, as each molecule is fragmented, the topomeric alignment and associated field is generated on-the-fly for each fragment and compared to the topomerically aligned field of a query molecule fragment. While the full fragmentation table scheme employed with the virtual library of application Ser. No. 08/903,217 may be employed, experience with fragmentations has shown that for medicinal type molecules the following fragmentation rule (which is a subset of the more general fragmentation method) produces meaningful fragments:

    
    
        “Break the molecule at acyclic bonds either singly or in pairs to generate sets of either 2 or 3 fragments respectively where each fragment must contain greater than a user specified number of heavy atoms.”
    
    


Assuming a setting that every fragment must contain at least three heavy atoms, FIG. 1 shows an example of how the rule is applied in a typical molecule (either a query molecule or a database molecule) to generate fragments. To generate the fragments, the whole structure is evaluated for each new fragmentation position. The two-piece fragmentations which will be performed are indicted by the thick lines. The two-piece fragmentations that will not be performed (because one of the resulting fragments contains less than three heavy atoms) are indicated by the thin lines. In this example, if, instead of requiring three heavy atoms, the user required five heavy atoms, then only the fragmentation between the two rings would be performed.
An example of a three piece fragmentation is shown in FIG. 2. Assuming again a setting that every fragment must contain at least three heavy atoms, the heavy lines indicate by arrows the two position in which the molecule would be fragmented into 3 fragments. The light lines indicate by arrows some of the three piece fragmentations that will not be performed because at least one of the fragments has fewer than three heavy atoms. If, instead of requiring three heavy atoms, the user required five heavy atoms, then no three-piece fragmentations would be performed.
At the present time, it has been found that generating three fragments is necessary when a two fragment scheme does not yield significant results. The three fragment scheme seems to find similar shapes that are sometimes missed in two fragment analysis. However, due to the higher computational overhead of three fragment searching, searches are first performed at the two fragment level. Four fragment searches may be necessary for some types of molecules, but at the time of filing the present disclosure, such situations have not been identified. Clearly the searching method of the present invention is not limited to the number of fragments which are generated but is generally applicable to as many fragments as the user wishes to consider.
Topomeric 3D Searching:
When analyzing molecules for shape similarity, it should be recognized that not all the elements of a molecule's shape may be required for proper interaction with a larger biomolecule. Perhaps in some instances, the entire shape is critical to the match. In other instances, only part of the molecule's shape may be critical to the match and other parts relatively unimportant. When comparing shapes of query molecules to those found in a heterogenous database, it is important to be able to compare not only the overall shape of the molecules, but also subparts. The method and software of the present invention permit many types of shape comparisons as will be discussed below.
Different heterogenous databases of compounds store compound structures in different formats such as SMILES, SLN, or an MDL format. Many software programs are available for interconverting the structures from one format to another. For the present application, the inventors use UNITY to convert compound information to SLN (Sybyl Line Notation) format. Compound information is then transferred to the CONCORD software program. CONCORD generates the three dimensional structure of the molecule. The starting point for topomeric searching of compounds listed in a heterogenous database are the CONCORD generated three dimensional structures of the database molecules and the query molecule. These structures are provided as input to the software programs set forth in the Appendices to the present disclosure.
The user specified fragmentation pattern (2 or 3 fragments and the number of included heavy atoms) is applied to the query molecule and the first database specified molecule. After each set of shape comparisons, the next database specified molecule is taken up in order. After the fragmentation patterns have been identified for each molecule (query or database), each fragment is aligned according to the topomeric rules.
In the preferred embodiment, the fragment is translated and placed into the grid so that the atom from which the “broken” acyclic bond extends into the fragment of interest is placed at the 0,0,0 coordinate. The “broken” bond (the attachment bond) is then directed along the X axis (standard topomer alignment) and the part of the molecule which is considered the fragment is aligned topomerically in the grided space. Alternatively, the atom in the fragment of interest which is connected to the acyclic bond which is “broken” is placed at the 0,0,0, position. This results in virtually insignificant differences in the topomer distances which are calculated.
Another feature of the present method is that a variable size grid region is used. Since some fragments are small and others large, the same volume of three dimensional grid space is not required to contain each fragment. Nothing is gained by placing a small fragment in a large grid space and only results in calculating an unnecessary number of extra grid location interactions. For the query molecule, the grid is adjusted to encompass the smallest region in which all the query fragments will fit. For database molecule fragments, the initial database molecule grid is one unit larger in all dimensions that the grid determined for the query fragments. The grid size is expanded by one unit in each dimension until the accumulated sum of the grid intersection points (starting with the query grid size and adding all the intersection points contained in each expanded grid) is greater than 10,000 or the grid has been expanded from its initial size by 11 units in each dimension. This procedure is followed since most computers, even those configured for molecular modeling, have a memory capacity which can be exceeded by allowing for unlimited grid size and number of intersection points. The grid size limitations are not required by the inherent method of the invention. Compression of the data from the thousands of data points in a large grid also aids in reducing the memory requirement for large grids. When a situation is encountered where the database molecular fragment extends outside of the maximum grid size, an “outside of the grid” factor is applied my multiplying the number of atoms outside the grid by the maximum interaction energy possible (typically 900) and adding that value as additional term in the root sum of squares similarity calculation. The use of dynamic grid sizing increases the throughput performance of the method considerably.
Whole Molecule Two Piece Comparisons:
As noted, for a two piece comparison both the query molecule and the database molecule are always split into just two pieces at each acyclic bond starting with the whole molecule each time. If there are 4 acyclic bonds and the heavy atom count matches the user selected value (default is typically=4), four two fragment pairs will be generated. As an example of the shape comparison, consider a query molecule which can only be broken at one acyclic bond to form fragments A and B. Consider also that a database molecule can only be broken at one acyclic bond into fragments C and D. Among the four fragments, there are two sets of comparisons possible: A:C & B:D, and A:D & B:C. A first comparison is made between: A:C and B:D. [In the actual calculation the squared differences in the field values between each grid location in each fragment are kept and the square root is only taken at the end of the comparison process.] Thus for the A:C & B:D comparison, a distance is determined as:

√{square root over ((A:C)2+(B:D)2)}{square root over ((A:C)2+(B:D)2)}

This value is retained for comparison. For the A:D & B:C comparison, a distance is determined as:

√{square root over ((A:D)2+(B:C)2)}{square root over ((A:D)2+(B:C)2)}

This value is compared to the value determined for the first A:C & B:D set and the lower value (greater similarity) retained. Thus, there are two comparison for each pair of molecules. It has been found that generally one will be significantly more similar than the other. The lower (more similar) value is retained and compared to the values obtained for the query against every other molecule in the database. Ultimately, the molecules in the database which are most similarly shaped to the query molecule will be determined by those with the smallest field difference.

As a further example consider a query molecule which can be broken at four acyclic bonds to form four two fragment pairs and a database molecule which can be broken at five acyclic bonds to form five two fragment pairs. this may be represented as:














Query
Database









A
I



B
J



C
K



D
L



E
M



F
N



G
O



H
P




Q




R










The first comparison will be A:I & B:J and A:J & B:I. A second comparison will be A:K & B:L and A:L & B:K. Similar comparisons will be obtained between each query fragment pair and each database molecule fragment pair. Of all the comparisons, the one having the smallest difference in field value will be kept for further comparison to the values obtained for all the molecules in the database. These comparison are whole molecule comparison because each fragment of the query molecule is compared to each fragment of every database molecule in sets of two (representing a complete molecule).
Whole Molecule Three Piece Comparisons:
If a three piece fragmentation scheme is employed the same shape comparison principles apply but are further complicated by the presence of the central fragment. In two piece fragmentation, each fragment has only one attachment bond which may be placed at the 0,0,0, grid coordinate. There is, therefore, only one topomeric alignment for the fragment. However, the central fragment in a three piece fragmentation will have two attachment bonds one each at the points were the two side fragments have been severed. There will, therefore, be two starting points for the topomeric alignment which will result in a different topomer shape of the aligned fragment. Each of these shapes must be included in the comparison.
As an example consider a query and a database molecule each which may be broken into three three piece fragmentations:














Query
Database









A
J



B
K



B′
K′



C
L



D
M



E
N



E′
N′



F
O



G
P



H
Q



H′
Q′



I
R










The primed fragments represent the second orientation of the central fragment of the three. Fields are calculated for all fragments as before. Considering just the first fragment set from both the query and database molecules the first set of distance comparisons are: A:J & B:K & B′:K′ & C:L and the distances is:

√{square root over ((A:J)2+(C:L)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:J)2+(C:L)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:J)2+(C:L)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:J)2+(C:L)2+[(B:K)2+(B′:K′)2]/2)}

The last term takes the average contribution of the center piece. Similarly, the other possible comparisons are calculated as:

√{square root over ((A:L)2+(C:J)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:L)2+(C:J)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:L)2+(C:J)2+[(B:K)2+(B′:K′)2]/2)}{square root over ((A:L)2+(C:J)2+[(B:K)2+(B′:K′)2]/2)}

From the two sets of comparisons, the one with the lower field difference (more similar) is retained for comparison. All the other comparisons between each three fragment set of the query and each three fragment set of the database molecule are calculated and the one with the lowest field difference is retained for comparison with those generated for all the other database molecules.

One further complication which arises with three piece fragmentation is that it is sometimes necessary to apply an attachment bond penalty to the calculated distance to reflect differences in the structure. Since there are two attachment bond points, the spatial relationship between those points will influence the shape of the whole molecule. However, considering just the fragments will not totally reflect the shape characteristics specified by the spatial relationship of the attachment points. This is an attempt to preserve the three dimensional structure of the whole molecule. A penalty value is thus added to the shape differences (increasing the apparent difference or similarity) to compensate. The penalty value is calculated as:

√{square root over ([(B:K)2+(B′:K′)2]/2)}{square root over ([(B:K)2+(B′:K′)2]/2)}

This penalty value is multiplied by an arbitrary factor depending on the user's belief in the significance of the structural difference, The penalty is initially set at 10 in the code but might be set as high as 100. For instance, as an example consider the ortho, meta, and para positional attachment bonds on a ring. The overall molecular shape will vary significantly if two side chains are in the ortho versus the para position with respect to each other. Accordingly, for the 1 atom difference of an ortho relationship, a penalty of 10 would be applied; for the 2 atom difference of a meta relationship, a 20 unit penalty would be applied; and for the 3 atom difference of a para relationship, a penalty of 30 would be applied. The point is that in determining the shape comparisons, a substituent can not just be moved around the ring and have it match without some penalty to reflect the difference in position.

For large molecules small changes in the number of atoms in the molecule is less likely to effect the overall shape than for small molecules. For effective shape comparisons, large structures need to be less sensitive to steric difference while small structures need to be more sensitive to steric differences. Experience has shown that there is a pivot point around 25 heavy atoms with structures considered large with more than 25 heavy atoms. Increasing the weighting of the steric contributions for small structures and decreasing it for larger structures has been found with experimental data sets to cut the number of false positives in half for small structures and allow more hits for large structures without eliminating many small structure hits.
Accordingly, for structures having more than 25 heavy atoms the steric field values calculated for each point in the grid may be decreased by as much as 33% (field values multiplied by 0.67). For structures having fewer than 25 heavy atoms the steric field values calculated for each point in the grid may be increased by as much as 100% (field values multiplied by 2.0). A non-linear multiple seems to work best.
In addition to using a variable grid size, another observation leads to a method of increasing the effectiveness and throughput of the searching methodology. It has been observed that for molecules which have a size difference of over +/−12 heavy atoms, there is little likelihood of finding molecules which match in shape. Consider a query with 20 heavy atoms and a database molecule with 33 heavy atoms. Since to start with there will be 13 atoms in the database molecule which will not be matched in the query, a large distance (dissimilarity) will already be found due to the missing atoms. The likelihood that all of the remaining atoms will lie in equivalent positions so that only the missing atoms will contribute to the difference in field values (and hence in similarity) is vanishingly small. Experimental runs on known data sets bears out this observation. Before any fragmentation is done, the difference in heavy atom size of the query and database compound is determined, and, if the difference is greater than 12 heavy atoms, the comparison is skipped.
Subset Searching:
As noted above, only part of the shape characteristic of many molecules may be responsible for the binding of those molecules to larger biomolecules. Accordingly, a search is desired which would find whether any part of the query molecule has the same shape as any part of the database molecule. This can be thought of as a partial fragment match. The method of this invention directly permits this type of search to be conducted. The query molecule is fragmented into two parts and the database molecule is fragmented into three parts in as many different ways as possible. For each possible three piece fragmentation you get:














Query
Database









E
A



F
B




C










In order to determine whether any part of the database molecule matches any part of the query the following comparisons are done:

















E:A
E:B



F:B
F:C



F:A
F:B



E:B
E:C










Since you are interested in locating any part of the database molecule which is closely similar in shape to all parts of the query molecule, the difference in heavy atom count exclusion which is applied to whole molecule searching is modified for subset matching. Instead of excluding the search if there is a +/−12 heavy atom difference, for subset searching the exclusion is not applied unless there is a +/−30 heavy atom difference.

Core Searching:
In some instances it is desirable to find another core of similar shape to a known core upon which a series of molecules may be built. For instance, suppose a patented series of compounds can be recognized as built upon a particular core. If that core can be replaced with a similarly shaped but chemically different core, it may be possible to construct an entirely new series of compounds active at the same site without infringing the patented series. To conduct this type of search the core and its two attachment bonds needs to be specified. How the searcher decides on the core structure is up to the searcher. The core is aligned in its two possible topomeric orientations and the fields calculated. The topomerically aligned field of only the central fragment of all possible three piece fragmentations of the database molecules are compared to the core fields as A:C & A′:C:














Query
Database









A
B



A′
C




D










Again, as before in the case of three fragment searching which involves a central fragment with two attachment positions, attachment penalties can be assigned to better characterize/distinguish the overall molecular shape based on where the attachment bonds are placed with respect to each other on the query core structure. For core searching, the penalty multiplier is typically set at 50. The molecules identified in the database which have central fragments generating the smallest values (greatest similarity) in the comparison to the specified core would be examined for possible use as cores.
Features:
As noted earlier, there may be some circumstances where the electrostatic field may be used in addition to the steric field to characterize the shape of a topomerically aligned fragment. A much more useful characterization has been implemented which extends ideas from pharmacophore modeling for use in searching heterogenous databases of compounds. It is well recognized that certain characteristic interactions of molecules in addition to shape play an important role in determining whether that molecule will bind to a larger biomolecule. Complimentarity of shape permits the molecules to approach each other closely enough for these interactions to take place. In pharmacophore modeling the presence and location of feature classes containing molecular characteristics thought important to the binding of the molecule is tracked as well as the distances and directions between the features. An absence of any given feature in a molecule or a different location is considered to significantly reduce the likelihood of that molecule's binding and, thus, typical pharmacophore modeling is an all or nothing proposition. Clearly, in the present methodology due to the topomeric alignment of fragments all distance and direction attributes of features present in the fragments are lost.
However, an alternative approach to incorporating the characteristic interactions in conjunction with the shape similarity matching described above has proven to generate an exceedingly powerful and accurate discovery methodology. The classic five feature classes are employed: positive charge, negative charge, hydrogen-bond-donating, hydrogen-bond-accepting, and aromatic. When present in either the query molecule or the database molecule, the features are assigned X,Y,Z point locations in the topomer alignment either centered on the relevant atom, or, in the case of aromatic rings, the centroid of the ring is specified. Generating the topomer conformation of a molecular fragment not only fixes the steric shape of that fragment, but is also fixes the Cartesian coordinates of each pharmacophoric feature contained within the fragment. The search strategy can be summarized as finding all the database molecule fragments which have features, similarly located in topomer space and similar in any other detailed feature property, that match each of the features in the topomerized fragments of the query structure.
In keeping with the distance definitions used for steric shape similarity, differences in features are defined with the same dimensionality as shape so that both shape and features can be used to characterize a fragment for searching. Feature by feature differences are also combined in a root sum square rather than a straight sum fashion. Thus, a second feature mismatch would not be as costly as the first one. To determine the feature “distance”, each of the pharmacophoric features in the query structure is considered in turn, by identifying the closest feature of the same pharmacophoric class in the database molecule fragment. If there is no such feature or if the nearest such feature is more than 1.5 Å distant, the dissimilarity sum of squares is increased by a maximum of 100×100 units. (Units are chosen to be commensurate with the steric shape units of kcal/mole-Angstrom3.) If there is a matching feature within 0.5 Å, the dissimilarity is set to zero. For a feature separation between 0.5 Å and 1.5 Å the dissimilarity penalty increment is obtained by linear interpolation between 0 and 100×100 unit values. Further, it is possible to scale/weight the feature contribution to increase or decrease its relative contribution with respect to the steric contribution to the observed similarity (distance).
Note that the use of the term “distance” with the feature searching methodology of the present invention is not meant to refer to an actual physical “distance” as considered in traditional pharmacophore techniques. For a two piece fragmentation the distance (similarity) between fragments is calculated as:














Query
Database









A
C



B
D










√{square root over ((A:C)FEATURES2+(A:C)STERIC2+(B:D)FEATURES2+(B:D)STERIC2)}{square root over ((A:C)FEATURES2+(A:C)STERIC2+(B:D)FEATURES2+(B:D)STERIC2)}{square root over ((A:C)FEATURES2+(A:C)STERIC2+(B:D)FEATURES2+(B:D)STERIC2)}{square root over ((A:C)FEATURES2+(A:C)STERIC2+(B:D)FEATURES2+(B:D)STERIC2)}

The cross terms for the A:D and B:C comparisons follow a similar definition as earlier. It has been observed that if the value of:

√{square root over ((A:C)FEATURES2+(B:D)FEATURES2)}{square root over ((A:C)FEATURES2+(B:D)FEATURES2)}

is too high, the distance will be large (little similarity) and the full calculation including the time consuming calculation of steric field can be skipped. This also increases the effectiveness and throughput of the method.

While the relative weight of each feature's contribution to the field can be varied, in the basic method, an attempt is made to match all features in a query with the nearest feature of the same class in the database molecule. This is similar to a pharmacophore type match, but there is no concern with matching interfeature distances in the topomeric conformation. Further, unlike standard pharmacophore searching, the user is able to assign adjustable penalties in the event that an exact match is not possible. For instance, a nearby spatial match of one type of feature might be more acceptable to the user than a nearby spatial match of another feature. The distance penalty for the spatially mismatched first feature could be set much lower than for a spatially mismatched of the second feature. The features method also permits handling of situations where a feature is present in a database molecule but not in the query molecule. In standard pharmacophore technique, this situation would lead to a total mismatch. However, in the present method the user can assign a distance (similarity) penalty for the absence of the match to the query, but need not totally ignore either the overall shape of the query or the contribution of the other features in judging the similarity of the structures.
Partial Feature Matching:
It is recognized that very frequently the binding of small molecules to receptors is highly dependant on the interaction between hydrogen-bond-donating and hydrogen-bond-accepting atoms. For partial feature matching, the search for charged groups and aromatic rings may be turned off. A large penalty (10,000 units) is applied for donors and acceptors which do not align. In addition, the number of donor or acceptor matches required can be varied. This capability is included since it is recognized that frequently only 2 or 3 groups are required to make a small molecule active. For partial feature matching, all the hydrogen-bond-donating and hydrogen-bond-accepting features are examined but only those generating the lowest 2 or 3 distances (including applicable penalties) across all (A:C, A:D, B:C, & B:D) the fragment comparisons for the compounds are used.
A further variation of the partial feature matching method considers the situation where the user determines that there is only one feature which is most important to match. If that feature is present and properly located, there is no penalty, the field differences are zero and the similarity is great. The flip side of single feature matching is that if the feature doesn't match a very large penalty is imposed to clearly yield a large difference (greater distance and low similarity).
Feature matching has been found to greatly increase the effectiveness of the heterogenous database searching since it compliments the shape specific searching. Use of both steric shape searching and feature searching of a topomerically aligned fragments has been found to be as good as or better than any equivalent 2D searching with fingerprints which has been, until now, the gold standard of searching technologies. In addition, the results of shape and feature similarity searching yields actual molecular structures which chemists recognize as being members of the same class of compounds. Also, unlike fragment searching, molecular structures are clearly identified which can serve as bases for continued development.
The method of the present invention for the first time permits the three dimensional searching of a heterogenous compound database for compounds that are likely to have the same biological activity as a query molecule. The results identify molecular structures having similar shape properties, and, when used with features, similar pharmacophoric properties. The identification of the structural fragments which contribute to the identified similarity provide an insight into the shape requirements of the receptor, and just as importantly, into likely additional molecular structures and corresponding shapes which will likely share the same activity. Thus, lead development is more straight forward from a knowledge of the relevant shape characteristics of the fragments provided by the method of this patent disclosure than from any two dimensional searching technique.
Output:
The most commonly used output reports the single best match between the query molecule and all molecules in the heterogenous database. The two or three piece fragment which was responsible for the match is also reported. A variation of the output, displays the fragment of the best hits and the query fragment that it matches. Once can also ask the system to list all hits with field differences less than some value; in other words a list of the most similar molecules.
The software code written in the C language contained in the Appendices implements all the capacities of the present invention. The CT_TOP.C code provides all the calculation functionalities. DBTOP.C contains the command line interface, the user inputs, code to read the input structures, calls to the CT_TOP.C routines, and output interface. CT_TOP.H lists all the required data structures used. The code needs to be compiled by a standard C compiler before being run as is well understood in the art. All together, all code necessary to fully disclose an enabling embodiment of the invention in the computational chemistry environment specified earlier is set forth in the Appendices.
From the proceeding description of the construction, generation, and searching of a heterogeneous database of molecules, it should be clear that there are many variations which may be employed and, having taught how to generate and search one specific embodiment, all equivalent embodiments are considered within the scope of this disclosure.
While the preceding written description is provided as an aid in understanding, it should be understood that the source code listings appended to this application constitute a complete disclosure of the best mode currently known to the inventors of the methods of heterogeneous database searching.
Thus, while this invention has been particularly described with reference to the drug lead identification art, it is clear that the validation of molecular structural descriptors and their use in selecting structurally diverse sets of chemical compounds can be applied anywhere a large number of compounds is encountered from which a representative subset is desired. Since the implications and advances in the art provided by the methods of this invention are still so new, the entire range of possible uses for the methods of this invention can not be fully described at the present time. However, such as yet identified uses are considered to fall under the teachings and claims of this invention if validated molecular structural descriptors are employed to characterize the diversity of molecules.

## Claims

1. A computer implemented method to search a heterogeneous compound database composed of molecules from different sources and syntheses, some known and some unknown, for molecules which are likely to have the same biological activity as a known query molecule comprising the following steps:
a) fragmenting a query molecule and database molecules according to a defined set of rules;
b) generating shape descriptors for the query molecule fragments and for the database molecules fragments; and
c) using the shape descriptors, comparing all combinations of query molecule fragments with database molecule fragments for each database molecule to identify the database molecule which has a shape similar to the query molecule; and
d) outputting the identity of the database molecule.

2. The method of claim 1 in which the output displays the fragment of the best hits and the query fragment that it matches.

3. A computer implemented method to search a heterogeneous compound database composed of molecules from different sources and syntheses, some known and some unknown, for molecules which are likely to have the same biological activity as a known query molecule comprising the following steps:
a) fragmenting a query molecule according to a defined set of rules;
b) topomerically aligning the query molecule fragments to generate a topomeric conformation;
c) generating the interaction energies between a probe and the atoms in the topomerically aligned query fragments at all intersection points in a three dimensional grid surrounding the aligned query fragments;
d) fragmenting a database molecule according to a defined set of rules;
e) topomerically aligning the database molecule fragments to generate a topomeric conformation;
f) generating the interaction energies between a probe and the atoms in the topomerically aligned database molecule fragments at all intersection points in a three dimensional grid surrounding the aligned database molecule fragments;
g) determining the similarity between query fragments and database molecule fragments by the root sum square differences in the field values; and
h) identify the molecule in the database most similar to the query molecule as that molecule having the smallest field value difference in its fragments; and
i) outputting the identity of the database molecule.

4. The method of claim 3 in which the output displays the fragment of the best hits and the query fragment that it matches.

5. A computer implemented method to search a heterogeneous compound database composed of molecules from different sources and syntheses, some known and some unknown, for molecules which are likely to have the same biological activity as a known query molecule comprising the following steps:
a) fragmenting a query molecule according to a defined set of rules;
b) topomerically aligning the query molecule fragments to generate a topomeric conformation;
c) generating the interaction energies between a probe and the atoms in the topomerically aligned query fragments at all intersection points in a three dimensional grid surrounding the aligned query fragments;
d) assigning features locations in the topomerically aligned query fragments;
e) fragmenting a database molecule according to a defined set of rules;
f) topomerically aligning the database molecule fragments to generate a topomeric conformation;
g) generating the interaction energies between a probe and the atoms in the topomerically aligned database molecule fragments at all intersection points in a three dimensional grid surrounding the aligned database molecule fragments;
h) assigning features locations in the topomerically aligned database molecule fragments;
i) determining the similarity between query fragments and database molecule fragments by the root sum square differences in the field values;
j) identifying all database molecule fragments which have features, similarly located in topomer space and similar in feature property, that match each feature in the query fragments; and
k) identifying the molecule in the database most similar to the query molecule as that molecule having the smallest field value difference in its fragments and smallest difference in feature; and
d) outputting the identity of the database molecule.

6. The method of claim 5 in which the feature contributions are weighted.

7. The method of claim 6 in which the output displays the fragment of the best hits and the query fragment that it matches.

8. The method of claim 5 in which only hydrogen-bond-donating and hydrogen-bond-accepting features are used.

9. The method of claim 8 in which the output displays the fragment of the best hits and the query fragment that it matches.

10. The method of claim 5 in which the output displays the fragment of the best hits and the query fragment that it matches.

11. A computer implemented method to search a heterogeneous compound database composed of molecules from different sources and syntheses, some known and some unknown, for molecular cores which are likely to have the same biological activity as a known query molecule core comprising the following steps:
a) specifying a known core and its two attachment bonds;
b) topomerically aligning the query core to generate a topomeric conformation;
c) generating the interaction energies between a probe and the atoms in the topomerically aligned query core at all intersection points in a three dimensional grid surrounding the aligned query core;
d) fragmenting database molecules into three fragments according to a defined set of rules;
e) topomerically aligning the central database molecule fragments generated by the fragmentation process of step (d) to generate a topomeric conformation;
f) generating the interaction energies between a probe and the atoms in the topomerically aligned central database molecule fragments at all intersection points in a three dimensional grid surrounding the aligned central fragments;
g) determining the similarity between query core and central database molecule fragments by the root sum square differences in the field values; and
h) identifying the core in the database most similar to the query molecule core as that core having the smallest field value difference molecule; and
d) outputting the identity of the core.

12. The method of claim 11 in which an attachment penalty multiplier is employed.

13. The method of claim 12 in which the output displays the fragment of the best hits and the query fragment that it matches.

14. The method of claim 11 in which the output displays the fragment of the best hits and the query fragment that it matches.