Statistical Modeling of Molecular Shape, Similarity, and Mechanism

Chapter 5

Statistical Modeling of Molecular Shape, Similarity, and Mechanism 1

Douglas R. Henry and A. Morrie Craig

2

1Molecular Design Limited, 2132 Farallon Drive, San Leandro, CA 94577 College of Veterinary Medicine, Oregon State University, Corvallis, OR 97331 Downloaded by GEORGETOWN UNIV on August 16, 2015 | http://pubs.acs.org Publication Date: November 14, 1989 | doi: 10.1021/bk-1989-0413.ch005

2

Current methods of molecular shape description and comparison have two major shortcomings: they are relatively nondirectional, and the superposition of molecules is usually a subjective procedure. A new method of shape description and comparison is described, which uses a binding moment, derived from fragment binding constants, to orient the molecules, followed by SIMCA modeling of the shape and similarity of the structures. This method is compared with other molecular shape descriptors of varying dimensionality. The shape descriptors are correlated with acute and chronic lung and liver toxicities of a set of pyrrolizidine alkaloids. The overall results, although not precise enough for quantitative prediction, demonstrate the value of using more specific and directional shape descriptors, coupled with a more rational method for molecule orientation and superposition. The interaction of a small molecule with an enzyme receptor, which is the basis for most therapeutic and many agricultural agents, proceeds in several stages (1). The first stage is a long-range through-solvent recognition and reorientation of the ligand, caused by interaction between the electrostatic fields of the ligand and the enzyme. In the next stage, the ligand approaches the receptor and becomes reoriented spatially and conformationally as it contacts and conforms to the surface of the active site. Finally, the ligand binds to the receptor, with possible perturbation of the receptor structure. This binding may be either a nonbonded interaction (reversible inhibition - the basis of action for many drugs), or it may involve formation of covalent bonds (irreversible and suicide inhibition - the basis for many agricultural agents). Depending on the stage of interaction, different aspects of a molecule's structure become important. At a distance, electrostatic field effects and the dipole moment of the structure are important; at first contact, steric effects play a role, and finally during the binding phase, more specific electronic and lipophilic properties come into play. Molecular shape plays an important part at all steps in the process. In solution, shape helps determine the dynamics of the molecule's movement and its interaction with solvent molecules. At first contact, shape is of paramount importance in limiting or enabling the interaction. In the final binding stage, the O097-6156/89W13-O070$06.00/0 © 1989 American Chemical Society

In Probing Bioactive Mechanisms; Magee, P., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1989.

Downloaded by GEORGETOWN UNIV on August 16, 2015 | http://pubs.acs.org Publication Date: November 14, 1989 | doi: 10.1021/bk-1989-0413.ch005

5. H E N R Y & C R A I G

Modeling ofMolecular Shape, Similarity, and Mechanism

shape of one part of the molecule may stoically interfere with a crucial rotation necessary to being a functional group into proximity with a specific region of the receptor. A precise and relevant description of the 'shape' of a molecule is difficult to calculate. Although several descriptors of molecular shape have been proposed, most fall short of being ideal. One way of viewing shape is from the standpoint of the dimensionality of the information the descriptor encodes. Kier has described a number of shape descriptors which are purely topological in nature (2). These do not depend on the 3D or even the 2D structure of the molecule, and they might thus be considered zero-dimensional measures of shape. Simple whole-molecule measures of size and bulk, such as the volume, the surface area, or any single dimension of the structure, might be considered one-dimensional measures of shape. Two-dimensional measures of shape include the cross-sectional shadow descriptors of Rohrbaugh and Jurs (3). Spherical harmonics provide a parameterization of three-dimensional shape which is especially useful for protein structure description (4). Overlap volume comparisons are the basis for the Molecular Shape Analysis (MSA) method of Hopfinger (5,6). This technique has been extended to include a quantification of the steric and electrostatic fields surrounding a molecule (7). A further refinement of field analysis, which merges statistical and molecular modeling techniques, is the COMparative Molecular Field Analysis method (COMFA) of Cramer (8). These latter approaches seek to encode information about more than just steric bulk or form. TTiey express multivariate information about the structure, so they might be considered multidimensional shape descriptors. A couple of problems exist with most of the shape descriptors mentioned above. The first is their nonspecific and nondirectional nature. In a single value, one cannot simultaneously encode information about both the direction and the form of the space the molecule occupies. Thus, a parameterized or at least a multivariate, description of shape is necessary. This is a feature of some of the descriptors mentioned, but it is usually not easy to interpret. A second problem deals with the orientation and superposition of molecules for shape comparison. Although methods have been published to automatically generate the best geometric fit to a set of matched points (9), in shape analysis this is almost always done subjectively, and it works best for structures which share a common parent substructure of atoms which can be superimposed in an obvious manner. There are many examples showing that chemical intuition can be misleading when orienting structures for comparison, especially if only 2D structures arc considered. A classic example is the flip of the pteridine ring between the binding of dihydrofolate and the binding of methotrexate to DHFR (10). To describe and compare molecular shape, and to correlate shape with biological activity, we sought a methodology that would not be limited by orientation and directionality shortcomings. We also wanted the method to be be implementable with existing modeling and data analysis software. The compounds we selected for analysis were a representative sample of toxic pyrrolizidine alkaloids, which are derivatives of the structures in Figure 1. These compounds produce lung and liver toxicity in horses and cattle which ingest the tansy ragwort (Senecio jacobaea). This problem has a multimillion dollar economic impact on agriculture in the Pacific Northwest, where the plant is abundant. Methodology

The structures and biological data were taken from an article by Culvenor and coworkers (11). There, 62 pyrrolizidine compounds were tested for acute and chronic lung and liver toxicity in young rats. The dose consisted of a single i.p.


72

PROBING BIOACTIVE MECHANISMS

HO

OH

OH

OH

RETRONECINE ^ O H

HELIOTRIDINE


HO DIOXO DERIVS.

DEHYDRO DERIVS.

Figure 1. Toxic pyrrolozidine structures. injection of the hydrochloride salt of the compound in neutral aqueous solution. Eight dose levels were used (0.025 mmol/kg to 3.2 mmol/kg). After four weeks, the rats were sacrificed, and microscopic examination was used to determine the level oftissuedamage. The authors defined five types of biological response: 1) peracute morbidity (death in 1 day or less); 2) acute morbidity (death in 1-7 days); 3) acute liver toxicity (centrilobular necrosis); 4) chronic liver toxicity (parenchymal megalocytosis); and 5) chronic lung lesions (intravascular and interstitial accumulation of mononuclear cells). The biological response was recorded as the lowest dose level producing the given lesions. A complication was that each of the compounds was nottestedat all dosage levels. The structures were quite diverse, ranging from those with simple alkyl ether attachments, to complex cyclized derivatives. The entire set could not be modeled and analyzed, so to obtain a representative sampling of the compounds, several molecular connectivity descriptors (12) were computed using the ADAPT program (13). These relate to the degree of branching in die structures, and they provide a convenient means of describing 2D structural diversity. The connectivity descriptors were combined with the biological data and a hierarchical cluster analysis was performed on range-scaled values (14). After analysis of the resulting dendrogram, 21 of the structures were selected (Figure 2). For the 21 structures, lung toxicity showed a fairly high correlation with acute death and liver toxicity. Otherwise, the biodata are fairly uncorrelated (Table I; P=peracute, A=acute, C=chronic). The 21 structures were modeled using the CHEMLAB-II system of Hopfinger and Pearlstein (15) and Allinger's MMP2 program (16). A rigorous conformational analysis was not performed for each structure; several low-energy conformers were identified for representative members of the set, and in the case of the macrocycles, their structures were compared with crystal structures for jacobine Table I. Biodata Correlation Matrix Death(P) Death(A) Liver(A) Liver(C) Lung

0.106 -.228 -.119 -.022

Death(A)

Liver(A)

Liver(C)

0.385 0.294 0.840

0.360 0.784

0.790



HEUOSUPINE

HEUOTRINE

AC-HEUOTRINE

LASIOCARPINE

RINDERINE


74


(17) and retrorsine (18). The conformer matching most closely was used in the analysis, and for building related structures. Standard molecular shape descriptors which were calculated included the following: Kier's K index, the Jurs shadow descriptors (6 descriptors) and length/breadth descriptors (2 descriptors), and Hopfinger's MSA volume descriptors (3 descriptors, using each structure in turn as a reference). The default conditions of die descriptor generating routines were used in each case. For comparison, some simple correlations were obtained between biological activity and the octanol/water partition coefficients, as calculated by the CLOGP program of Leo (19). An example is seen in Figure 3, which shows an apparent parabolic relationship between the log(l/minimum dose for acute death) and log P. The two outlier structures may have poorly calculated log P values (the epoxide ring and the conjugated system are special features), or they may undergo special metabolic transformation. A Binding Moment Approach to Molecular Orientation. Except for the Kier index, each of the conventional molecular shape descriptors requires a standard orientation of the structure. In the case of the shadow and length/breadth descriptors, this is along principal axes of the molecule. In the case of the MSA descriptors, the structures were aligned by matching the following three atoms:


2

As an alternative to specific atom matching, a moment of binding was calculated using the fragment binding constants of Andrews (20). The binding moment was computed in a manner similar to the calculation of a dipole moment from partial residual atomic charges. The partial charge values for the nonhydrogen atoms in the structure were replaced by binding constant values. In the case of carbonyl groups, the binding constant value was positioned on the oxygen atom. Alternatives would have been to split the binding constant value, or to generate a pseudoatom at the centroid of die atoms in the fragment. The calculated moment was scaled to unit length and stored as pseudoatoms with the structure. Figure 4 shows a plot of the MMP2-modeled structures of senecionine and dehydrosenecionine, with circles representing the fragment binding constants, and the binding moment displayed as a vector. Also shown are circular profiles of the biological activity. The length of any given spoke in the profile is proportional to the level of activity, and here one sees that the compounds have quite different spectra of biological activity. To superimpose structures for comparison, the corresponding ends of the binding moment vector, plus the ring nitrogen atom, were used asreferencepoints. In many cases, this gave superpositions similar to those obtained using the nitrogen and the two oxygen atoms, in the MSA analysis. In some cases, as shown in Figure 5, there was a reorientation of the structure being superimposed, sometimes amounting to an almost 90 degree rotation about the vertical axis. It was interesting to note that these differences in orientation often occurred for structures which had quite different toxicity spectra. Of course, other moments, such as lipophilicity and dipole moments, could be matched for superposition as well, though this was not done in this case. SIMCA Modeling of Molecular Shane. To model the shapes of the molecules in a manner that would express direction as well as bulk, and allow comparison of


5.

HENRY & CRAIG

Modeling ofMolecular Shape, Similarity, and Mechanism L o g P V e r s u s A c u t e Death


3. -\

C a l c u l a t e d L o g P (Jurs)

Figure 3. Correlation of calculated log P (octanol/water) values with acute morbidity of pyrrolizidine alkaloids. Structures of two outlying compounds are shown.

X Y Z Rotations: 0.0

0.0

0.0

X Y Z Rotations: 0.0

0.0

0.0

Figure 4. Modeled structures of senecionine and dehydrosenecionine, showing circles for fragment binding constants and the binding moment vector. Circular profiles of biological activity are also shown. In Probing Bioactive Mechanisms; Magee, P., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1989.


76


different structures, we chose the SIMCA (Soft Independent Modeling by Class Analogy) technique (21). This method generates disjoint principal component models of clustered points in multidimensional space. With SIMCA, it is possible to describe and parameterize spherical (zero components), linear (one component), planar (two components), and box-shaped regions (three and higher components) of space. The usual limitations apply in terms of the statistical degrees of freedom. In this case, the points were the nonhydrogen atoms in the structure. The variables were simply the Cartesian coordinates of the atoms. It is also possible to include other atomic properties, such as electrostatic and lipophilic characteristics, if necessary, to generate an extended definition of molecular shape. This was not done here, since we were primarily interested in conventional shape description. SIMCA is not a clustering techniqe. Thus, it was first necessary to cluster the atoms in each structure. This was done using hierarchical clustering in Cartesian space, with ADAPT. A sample cluster dendrogram for senecionine is shown in Figure 6. The number of clusters for each structure was manually determined, though we could have applied one of many cluster validation techniques (22). Between three and five clusters were chosen for each structure. As one would expect for bonded atoms, the clustering followed bonding patterns, since bonded atoms are closest together. Once die clusters were selected, SIMCA analysis was performed on unsealed coordinate data, treating each cluster as a separate class of points. Default cross-validation techniques led to between one and two principal component models for each cluster, A schematic representation of the SIMCA models for the atoms in senecionine is seen in Figure 7. Comparison of Structures. To compare the structures and generate shape descriptors, the reference compound approach was used. Each compound was used in turn as a reference structure. Every other structure was superimposed on the reference using the binding moment technique described previously. Then, each nonhydrogen atom in the sur^rimposed structure was treated as a new "observation" in the analysis. The atom was fit to the principal component models of each cluster of atoms in the reference structure. The atom was 'assigned* to the cluster itfitbest, and both the standard deviation value and the F-value for the fit was recorded. The final descriptor value was taken as the average standard deviation or the average F-value for all the atoms placed in a given cluster. This technique incorporates the statistical variation of a given model, so that a structure will not 'fit' its own models perfectly. The fit of a structure to its models, however, defines a lower range of values for comparison with other structures. As values increase beyond this lower range, it indicates a poorer shape comparison. The fact that disjoint models are used means that directionality can be expressed. When senecionine is used as the reference structure, one obtains the standard deviation and F-values shown in Table II. Values are shown for senecionine fitted to its own models, and for a similar structure (seneciphylline) and a dissimilar one (dehydrosenecionine). As this table shows, the F-values are largo and show a wider variation than the standard deviations. They generally correlated with the standard deviation values, but in some cases (cluster 2 in Table II), there were discrepancies. Because of their wider span, and because they generally accorded with graphical comparison of the structures, the F-values were used in the correlations with biological activity and for comparison with other shape descriptors. Results and Discussion In general, the various different shape descriptors in this study were not intercorrelated (R values 0.7 or less). Simple correlations of the SIMCA F-value descriptors with the various biological response variables yielded rather poor results.



5.

H E N R Y & CRAIG

Modeling of Molecular Shape, Similarity, and Mechanism 77

Figure 5. Superposition of senecionine and dehydrosenecionine obtained by matching binding moment vectors and ring nitrogen atoms. The dehydro compound is rotated almost 90 degrees relative to senecionine. T r a i n i n g Set: 0

i

—

Metric: Euclidean Fusion:

Furthest N e i g h b o r

Clustering of Patterns |

D A N CLASS

I

2 2 2 2 1 1 1 1 1 1 2 1 1 3 5 2 4 6 7 8 9 1 1 1 2 0 4 3 5 6 4 2 3 1 8 9 7 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

1 1 1 1 1

Figure 6. Dendrogram showing atom clusters for senecionine. Three clusters were selected to represent the structure. In Probing Bioactive Mechanisms; Magee, P., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1989.



w

en

S

î

η

§

ο

δ

sΟ

00

5.

HENRY & CRAIG

Table n.

Modeling of Molecular Shape, Similarity, and Mechanism Average and Maximum Standard Deviation and F-value Descriptors with Senecionine as Reference (maximum in parentheses) Fitted Compound

Senecionine

Seneciphylline

Dehydrosenecionine

1 (sd) 1 (F)

0.218 (0.510) 0.469 (1.921)

0.643 (1.313) 3.590 (13.43)

0.873(1.799) 5.020 (20.45)

2 (sd) 2 (F)

0.557 (1.390) 0.625 (2.350)

0.907 (1.610) 1.680(4.030)

0.505 (1.322) 5.570 (26.24)

3 (sd) 3 (F)

0.230 (0.530) 0.547(1.928)

0.663 (1.245) 1.500(3.501)

0.958 (2.458) 4.436(11.56)


Cluster

Correlation coefficients ranged from 0J6 or lower to 0.7, depending on the reference compound that was used. With only a few exceptions, similar results were obtained for die other descriptors. Since the biodata and the shape descriptors in this study are both multivariate, canonical correlation analysis was selected to provide a single overall measure of correlation between molecular shape and biological activity, for comparison of the various shape descriptors. In canonical correlation analysis, the combination of the predictor variables is found, which correlates highest with any possible combination of the response variables (23). A similar approach is taken in Partial Least Squares (PLS) analysis (24). Table III shows, for the 21 structures in the analysis, the first canonical correlation coefficients relating the shape descriptors with the five measures of biological activity. In the case of a single descriptor (Kier's K index, for example) the canonical correlation coefficient is the same as the simple correlation coefficient, so univariate and multivariate correlations can be compared directly. The correlations in this table are the highest that were observed in the analysis. The activities in the last column are the biological response variables that were most highly associated with thefirstcanonical variable of the biological data. The best correlation was observed with Jurs* shadow descriptors. This may be partly a result of the larger number of descriptors used, since the simplest topological and one-dimensional descriptors show the lowest levels of correlation. The MSA and SIMCA descriptors are directly comparable. The SIMCA descriptors have somewhat higher correlations with activity, though the differences may not be statistically significant. They are computed with about the same amount of effort as the MSA descriptors; in addition, the SIMCA descriptors have directionality, which could allow a researcher to determine which part of the molecule is responsible for the shape differences and presumably, the differences in activity. The quality of the biodata was not really high enough to associate particular types of activity with certain regions of the molecules. In the future, we hope to be able to accomplish this with newer in vitro assays of pyrrolizidine toxicity. 2

Conclusions As the technology for computer simulation of 3D molecular structures improves, we are able to generate ever more accurate models of structures and activity. A recurrent trend in QSAR science has been the shift from the more general to the


80


Table HI.

Vars.

Canon. Corr.

Highest Correlating Activity

Kier 2 K

1

0.50

Liver(C)

Jur's Shadow

6

0.86

Death (P)

JUTS'

2

0.42

Death (P)

MSA Senecionine Usaramine Madurensine Average

3 3 3

0.72 0.74 0.64 0.66

Liver (A) Liver (C) Lung

SIMCA Retrorsine Rinderine Madurensine Average

3 3 4

0.80 0.75 0.78 0.72

Liver (A) Liver (C) Lung

Descriptor Set


Canonical Correlation Results

L/B

more specific (as in proceeding from bulk log P to % values, and in moving from topological to geometrical descriptors). More accurate models of structure require more specific shape descriptors to adequately encode relevant information for correlation with biological activity. As we have shown, there are advantages to multivariate descriptions of molecular shape based on statistical modeling of the structure. Such models are capable of expressing both the amount and the direction of shape differences. We are presently investigating the use of lipophilic and electrostatic extensions of our definition of molecular shape. Literature Cited 1.

Franke, R. Theoretical Drug Design Methods; Elsevier: New York, 1989; pp 316-322. 2. Kier, Lemont B.Quant.Struct-Act.Relat. 1986, 5, 11-12. 3. Rohrbaugh, R. H.; Jurs, P. C. Anal. Chem. 1987, 59, 1048-1054. 4. Max, Nelson; Getzoff, Elizabeth D. IEEE Comput. Graph. Appl. July 1988; pp 42-50. 5. Hopfinger, A. J. J. Am. Chem. Soc. 1980, 120, 7196-7206. 6. Hopfinger, A. J.; Burke, Benjamin J. In QSAR: Quantitative Structure-Activity Relationships in Drug Design; Fauchere, J. L. Ed.; Alan R. Liss: New York, 1989; pp 151-159. 7. Hopfinger, A. J. J. Med. Chem. 1983, 26, 990-996. 8. Cramer, R. D.; Patterson, D. E.; Bunce, J. D. J. Am. Chem. Soc. 1988, 110, 5959-5967. 9. Danziger, D. J.; Dean, P. M. J. Theor. Biol. 1985, 116, 215-224. 10. Kuyper, Lee F. In Computer-Aided Drug Design; Perun, Thomas J; Propst, C. L. Eds.; Marcel Dekker: New York, 1989; pp 337-338. 11. Culvenor, C. C. J.; Edgar, J. A.; Jago, M. V.; Outteridge, A.; Peterson, J. E.; Smith, L. W. Chem.-Biol. Interactions 1976, 12, 299-324. 12. Kier, Lemont B.; Hall, Lowell H. Molecular Connectivity in Drug Research; Academic: New York, 1976. In Probing Bioactive Mechanisms; Magee, P., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1989.

5.

H E N R Y & CRAIG

13. 14. 15. 16.


17. 18. 19. 20. 21. 22. 23. 24.

Modeling of Molecular Shape, Similarity, and Mechanism 81

Stuper, Andrew J.; Brugger, William E.; Jurs, Peter C. Computer Assisted Studies of Chemical Structure and Biological Function; Wiley: New York, 1979. Romesburg, H. C. Cluster Analysis for Researchers: Lifetime Learning Press: Redwood City, CA, 1984. Pearlstein, R. A. In Chemlab-II Reference Manual; Molecular Design Limited: San Leandro, CA, 1988. Sprague, Joseph T.; Tai, Julia C.; Yuh, Young; Allinger, Normal L. J. Comp. Chem 1987, 8, 581-603. Rohrer, D. C.; Karchesy, J.; Deinzer, M. Acta. Cryst. 1984, C40, 1449. Coleman, P. C.; Coucourakis, E. D.; Pretorious, J. A. S. Afr. J. Chem. 1980, 33, 116. Leo, A. In Medchem Software Manual - Release 3.52; Daylight Chemical Information Systems: Irvine, CA., 1987; Chapter 14. Andrews, P.R.; Craik,D. J.; Martin, J.L. J. Med. Chem. 1984, 27, 1648-1657. Wold, Svante Pattern Recognition 1976, 8, 127-139. Milligan, Glenn W.; Cooper, Martha C. Psychometrica 1985, 50, 159-179. Morrison, Donald F. Multivariate Statistical Methods; McGraw-Hill: New York, 1976; pp 259-263. Wold, Svante; Geladi, Paul; Ebensen, Kim; Ohman, Jerker J. Chemometrics 1987, 1, 41-56.

RECEIVED

August 2, 1989


Statistical Modeling of Molecular Shape, Similarity, and Mechanism

Recommend Documents