Chem. Res. Toxicol. 2006, 19, 1533-1539
1533
Top-Priority Fragment QSAR Approach in Predicting Pesticide Aquatic Toxicity Mose’ Casalegno,*,† Guido Sello,‡ and Emilio Benfenati† IRFMN, Mario Negri Institute for Pharmacological Research, Via Eritrea 62, 20157 Milano, Italy, and Dipartimento di Chimica Organica e Industriale, UniVersita’ degli Studi di Milano, Via Venezian 21, 20133 Milano, Italy ReceiVed August 7, 2006
In the framework of pesticide risk assessment, a fragment-based QSAR approach is presented to correlate LC50-96 h acute toxicity to the rainbow trout (Oncorhynchus mykiss). While there are other fragmentbased modeling routes, our approach exploits the possibility of prioritizing fragments’ contributions to toxicity. On the assumption that one fragment might be mainly responsible for the molecular toxicity, we developed a three-stage modeling strategy to select the most important moieties and to establish their priorities at a molecular level. This strategy was tested on a heterogeneous dataset containing 282 pesticides, collected under the EU-funded project Demetra. Quantitative toxicity prediction yielded good results for the training set (R2TR ) 0.85) and the test set (R2TS ) 0.75). The advantages and limitations of the current priority strategy are examined. Introduction Pesticides have become essential products in everyday life. Modern pesticides are used to protect plants and crops against pests, to control overgrowth of undesirable plant species, and to protect public health from disease vectors such as mosquitoes, ticks, rats, and disease-causing organisms. Significant social and economic benefits are related to the use of pesticides: fatal infections are kept under control, and safe food is available all year round. Pesticides contain a variety of organic and inorganic compounds. By nature, they are poisonous, and while they are useful against undesirable target organisms, widespread and mismanaged use can pose risks to nontarget organisms and the environment. For this reason, pesticide approval, registration, and marketing are subject to strict control by regulatory authorities both in the U.S. and EU (1, 2). In both legislative frameworks, risk assessment procedures have been implemented to measure the impact of pesticides on different environmental compartments (3-5). Pesticide risk assessment is a wide and complex process. The toxicological properties of each pesticide must be fully evaluated toward several nontarget animals and plants. All breakdown products from environmental degradation or metabolic transformations must also be evaluated. Inclusion of pesticide transformation products steeply increases the number of compounds to be tested and requires massive efforts in terms of money and testing animals. For this reason, the adoption of QSAR methods has been considered in both U.S. and EU official documents (6, 7) as a valid alternative to animal testing. QSAR can help identify compounds of potential concern during the preliminary phase of risk assessment. Those estimates serve as a complementary tool to support decision-making on whether further testing is necessary, to fill experimental data gaps, and to set testing * To whom correspondence should be addressed. Tel: ++39-0239014586. Fax: ++39-023546277. E-mail:
[email protected]. † Mario Negri Institute for Pharmacological Research. ‡ Universita’ degli Studi di Milano.
priorities. Implementation of QSAR strategies helps avoid unnecessary animal testing, boosting cost- and time-effectiveness in risk assessment. In this context, fragment-based QSAR models are of potential interest, in view of the possibility of identifying biologically active substructures responsible for the pesticidal activity. In contrast to holistic approaches, which involve descriptors referring to the whole molecular structure, fragment-based methods investigate structure-activity relationships at the substructural level. This means structure and toxic effects can be related explicitly in a transparent and intuitive way. Once biologically active moieties (toxicophores) have been extracted and collected, unknown molecules could be inspected, checking for their presence. Test candidate molecules might be new pesticides, yet to be manufactured, or metabolites or degradation products from transformation or degradation processes. The aim is always to identify potentially hazardous compounds in advance, following the assumption that substances containing the same toxicophore are likely to cause similar toxic effects (8). Then, if the active moiety remains intact during pesticide degradation, the transformation product may exert its toxic potency the same way as the parent compound. Conversely, if the toxicophore is removed or altered after degradation, it is unlikely that the transformation product will show the same activity (9). New active fragments can be added to the original pool as soon as new data become available, to include all the main pesticide classes, with a view to large-scale pesticide evaluation. In principle, this idea offers an opportunity to rationalize pesticide toxicity. In practice, however, its potential has been only exploited marginally. Among the QSAR models for correlating pesticide toxic effects toward nontarget species so far reported (10-17), only a few rely on fragment-based architectures (18-21). The lack of applications in this field is mainly due to the limited availability of adequate pesticide data sources. Fragmentbased approaches demand large, unbiased training sets to provide accurate predictions. A limited number of data points might lead
10.1021/tx0601814 CCC: $33.50 © 2006 American Chemical Society Published on Web 10/26/2006
1534 Chem. Res. Toxicol., Vol. 19, No. 11, 2006
to a reliable model with limited applicability (22, 23). However, increasing the data with the aim of enlarging the applicability domain does not necessarily guarantee accurate predictions. A number of statistically significant fragments must be collected to ensure objective investigation and reliable model building. Meeting this requirement is possible as far as homogeneous datasetssi.e., containing structurally similar compoundssare considered. However, in QSAR practice, datasets are often small and heterogeneous. This results in a multitude of fragments describing unique chemical environments, with the risk of incorrect identification of biologically active moieties. The quality of the biological data provided, with the chemical structures, is also important for the conclusions that will be drawn from the model’s results. The lack of data might pose an insurmountable barrier to the adoption of fragment-based methods. With the absence of adequate toxicological data, no efforts have gone into addressing some important methodological issues. In QSAR practice, fragments have often been treated as traditional holistic descriptors in model building. Most of the statistical approaches successfully applied to holistic descriptors, such as multilinear regression and neural networks, have been exploited in developing fragment-based models. Additive linear equations have become a natural choice among investigators, owing to their simplicity and ease of implementation. Here fragments behave as individual, noninteracting entities. The toxicity is estimated by summing the individual contribution of each fragment in the molecule. Fragment contributions, expressed as real coefficients, are optimized by minimizing a cost function, typically the error between predicted and observed responses. It follows that all fragments have the same chance of being used in estimating the toxicity. This route differs considerably from that suggested above for pesticides, where only toxicophores play a leading role in determining the toxicity. In this case, fragments lose their equality, active moieties gaining importance at the expenses of inactive ones. The interactions between fragments cannot be neglected, because the presence of two or more toxicophores in a molecule raises the issue of prioritizing their toxic effects. These considerations pose a question about the validity of additive linear equations in modeling pesticides and led us to develop a method where fragments’ contributions to the toxicity are not equal but prioritized. We assume that, out of a set of molecular fragments, one is mainly responsible for the molecular biological activity. According to this fundamental assumption, this work describes a three-stages modeling route to pesticides resembling that outlined above. The strategy was applied to 282 pesticides collected under the EU-funded project Demetra (24). The acute toxicity data (LC50-96 h) for rainbow trout (Oncorhynchus mykiss) were carefully selected by following a robust and standardized procedure (25). This gave us a unique opportunity to overcome the limitations arising from the use of poor-quality biological data. At the same time, this dataset provided a very challenging test bed for our method, as indicated by the number of pesticide classes included, more than 20. The whole dataset was partitioned in training and test sets, with respectively 240 and 42 molecules, using the partition selected in the Demetra project (24). Training molecules were used to extract the priority rules and to check the applicability of the priority scheme. Test molecules were submitted as an external set, to measure the model’s performance and check the validity of our fundamental assumption.
Casalegno et al.
Materials and Methods The first stage of our route involves extracting the most important fragments from the training set. All 282 molecules are broken down into small substructures, called atomic centered units (ACUs). A numerical criterion based on the training set toxicity data is then applied, to assign one fragment to each training set molecule, called top-priority fragment (TPF), made up of one or more ACUs. Once the TPFs are extracted, in the second stage, we evaluate their priority relationships. A “priority matrix” is used to retrieve all priority relationships. In the last stage of our modeling route, molecules are submitted and inspected to check for the presence of TPFs. If more TPFs are found, the information in the priority matrix is exploited with the aim of identifying the one(s) having highest priority. The quantitative prediction is based on average fragment toxicity. For the 282-pesticide toxicity data for the rainbow trout (Oncorhynchus mykiss) we considered the experimental acute toxicity LC50-96 h, that is the median lethal concentration after 96 h of exposure, expressed in mmol/L. For modeling purposes, the log(1/LC50-96 h) was used instead of the original LC50-96 h. 2D chemical structures were drawn for all compounds in MDL MOL format. In a preliminary stage, the MDL MOL files were preprocessed, so as to optimally restructure the chemical information. This procedure was described in a previous work (26), with the difference that, in the present case, the program was instructed to detect fivemembered heterocyclic aromatic systems, either single or condensed, in addition to six-membered ones. To all bonds belonging to a five- or six-membered aromatic ring, a bond order of 4 was assigned instead of the original value of 1 or 2. In addition, in condensed aromatic systems the bond order of common bonds was increased by 1. For example, a bond order of 5 would be assigned to the bond common to both naphthalene aromatic rings. It should be noted that ACUs in nonaromatic cycles have not a special code; for example, there exists only one ACU for the fragment CH2CH2-CH2 both in a cyclohexane and in a chain substructure. The molecular representation given above provided the input to generate ACUs. Atomic centered substructures have been described in the literature (27); these are simply made by a central, or parent, atom and its nearest neighbors. To generate them, molecules were broken down through the following procedure. For each molecule, one parent atom was picked up in turn, and then all its nearest neighbors were selected. The resulting substructure was then isolated from the rest of the molecule by cutting all bonds connecting it to other atoms. The following chemical features were used to define each ACU: (1) the elements making up the constituent atoms (parent one and its nearest neighbors), such as P, O, C, etc.; (2) the bond orders referring to all bonds joining the ACU constituent atoms; (3) the bond orders referring to all bonds joining nearest neighbors with other atoms. Application of the breakdown rule yielded a total number of ACUs equal to the number of atoms for each molecule. Figure 1 illustrates the result of this procedure for a small molecule (CAS No. 16752-77-5) containing 10 atoms. A total of 10 ACUs, one for each parent atom (bold typed), were obtained. All ACUs were collected and compared so as to eliminate redundancies and to compute fragments occurrences for all 282 pesticides. This gave Ntot ) 730 different ACUs, describing the entire dataset. An integer number was assigned to each ACU, from 1 to Ntot. Then, for each molecule, a binary fingerprint array, B, was filled, accounting for the jth fragment presence (B ) 1) or absence (B ) 0). At this point, the complete dataset was divided into two subsets, a training set (240 molecules) and a test set (42 molecules), giving three classes of ACUs: those occurring only in training molecules, those common to both training and test molecules, and those occurring only in test molecules. The ACUs present only in the test molecules were obviously discarded; in addition, we discarded
Top-Priority Fragment QSAR Approach
Chem. Res. Toxicol., Vol. 19, No. 11, 2006 1535
Figure 1. Example of break down procedure. ACUs parent atoms are bolded.
the ACUs that appeared only in one training compound because their use in prediction should be limited. This reduced the total amount of ACUs to Neff ) 334 effective units. Only training set molecules were considered in building the prioritization scheme, consisting of three sequential stages: (1) Preassignment. This preliminary step involved implementing a criterion to identify important fragments. Usually these are selected out of a pool of predefined fragments, but here we followed a different route. We adopted a numerical criterion to assign one top-priority fragment (TPF) to each training set molecule. To assign a toxicity to each ACU, we used training set toxicological data. For each ACU we selected all the compounds containing it; then we summed the observed molecular toxicities (Tobs) of the selected compounds and divided the sum by the number of molecules to compute the average toxicity (Tave) of the current ACU. This procedure gave a simple way to assign to each ACU an activity coefficient equal to its arithmetic mean toxicity in the training set. Then for each training set molecule, the absolute difference between the observed molecular toxicity and the average ACU toxicity was minimized, i.e. |Tiobs - Tjave|
for j ) 1, Neff
(1)
where Tiobs is the observed toxicity for the ith molecule and Tjave is the average toxicity for the jth ACU. In each molecule, all ACUs minimizing the residual expressed by residual (1) were grouped into a single fragment, the TPF for that molecule. At the end of the procedure, one TPF, comprising one or more ACUs, was thus assigned to each training set molecule. According to residual (1), all ACUs grouped within the same TPF had the same average toxicity. This permitted us to associate an average toxicity with each TPF equal to that of any of the constituent ACUs. Looking at (1), it is easy to see that the criterion assigned top priority to ACUs whose Tave were closest to the observed molecular toxicity. Within our main assumption, this was the best choice from the quantitative point of view. In other words, the best fragment for estimating the molecular toxicity is the one accounting for the largest part of the experimental value. The number of TPFs resulting from this procedure, hereafter indicated as Ntpf, corresponded to the toxicophores considered in the next stage. (2) Setting and Retrieving Reciprocal Priorities. The preassignment procedure provided two results: the selection of a pool of Ntpf toxicophores and the assignment of a toxicophore to each training set molecule. Both were done using known toxicity data. Now we need a general procedure for identifying the TPF when the toxicity is not known, i.e., in a test molecule. To identify the TPF of the test molecule we use the TPFs of the training set compounds. First, TPFs have to be retrieved within the molecule under investigation. Then, if only one of the known TPF is found in a test compound, we assign it to the molecule. Otherwise if two or more TPFs are present we select the one best suited to represent the molecule’s toxicity. To do this, we introduce
the concept of priority relationship among TPFs. This is a mathematical statement that serves to rank two TPFs. Pairwise priorities can be conveniently stored in a priority matrix, M, a square matrix, consisting of Ntpf × Ntpf elements, whose generic elements MJ,JJ retrieve the priority relationship between the fragments J and JJ. To fill the priority matrix, the roles of fragments J and JJ were investigated in the molecules where they were concurrently contained. Four different types of relationship were considered. Below Table 1 summarizes them and the conditions under which they were assessed. A capital letter (first column), a type (second column), and a conditional statement (last column) were the main attributes reported for each relationship. The conditional statement lists the conditions to be satisfied to assign a label character to each matrix element, MJ,JJ. An example might help illustrate the retrieval procedure. Suppose that J and JJ concurrently occur in six molecules. Now, suppose that fragment J has been assigned to five of them as TPF. If JJ has never been assigned as a TPF within this subset, J “wins” over JJ. The corresponding matrix element, MJ,JJ, is therefore “W”. At the same time, the matrix element MJJ,J is set at “L”, since the fragment JJ “loses” priority with respect to J. However, if JJ has top priority for one molecule, the relationship is “conflicting” because both fragments have top priority at least once. A conflicting relationship arises every time both fragments are found to have top priority, regardless of the number of molecules in which they play this role. According to the rules reported in Table 1, in this case, both matrix elements MJ,JJ and MJJ,J are set as “C”. The null, “V”, relationship occurs when fragments J or JJ have top priority in none of the six molecules. It might also happen that two fragments are never concurrently present in the six molecules. In this case the molecular subset containing both fragments is empty, and their mutual relationship is, therefore, “void”. The conditional statements reported above are mutually exclusive. It is not possible for two TPFs to have more than one kind of relationship. At the same time, the conditional statements account for all possible relationships between two fragments. An excerpt of a typical priority matrix after the retrieval procedure is illustrated in Figure 2. The numbers aligned along both axes served to identify the TPFs. The diagonal elements, having no role in the priority scheme, were labeled with a dash. For example, the element M1,2 states that TPF 1 wins over TPF 2. Its symmetrical counterpart, M2,1 is exactly the opposite because TPF 2 loses. (3) Postassignment (Prediction). Postassignment is the last stage of our modeling strategy, where the information retrieved in the priority matrix is exploited for predicting the molecular toxicity. The goal is to assign one TPF to a molecule and to predict its toxicity. The current molecule is initially inspected, looking for TPFs, as only TPFs play a role in predicting toxicity. Depending on the number of TPFs found, three different situations could be expected: (A) No TPFs Found. In this case, no prediction is possible for the molecule under investigation. The molecule contains no active moieties and, therefore, lies outside the domain of applicability of the model. (B) One TPF Found. As this is the sole fragment, it maintains highest priority. The average TPF toxicity is directly used to estimate molecular toxicity. (C) Many TPFs Found. In this case, the priority matrix comes into play. Priority relationships referring to these TPFs are extracted from the priority matrix. For each TPF, the winning, losing, and conflicting scores are computed, obtained by summing the number of winning, losing, and conflicting relationships, respectively. The TPF with the highest winning score is the fragment having highest priority. If two or more TPFs have the same winning score, the one with the lowest losing score has highest priority. Finally, if two or more TPFs also have the same losing score, the lowest conflicting score is used to assign top priority. At this point, no more selection criteria are available to further discriminate TPFs.
1536 Chem. Res. Toxicol., Vol. 19, No. 11, 2006
Casalegno et al. Table 1. Priority Relationships
MJ,JJ
relationship
W
winning
L
losing
C V
conflicting null, or void
conditional statement (within the subset containing both) J is the top-priority fragment in at least one molecule JJ is never the top-priority fragment JJ is top-priority fragment in at least one molecule J is never the top-priority fragment both J and JJ are top-priority fragments at least once J and JJ are never top-priority fragments, or they are never concurrently present Table 2. Four Molecules Not Correctly Reassigned
Figure 2. Excerpt of a priority matrix.
If one top-priority fragment is finally picked up, it is assigned to the molecule and its average toxicity is used as the predicted toxicity. If, however, the selection ends with more than one TPF having the same priority, all are assigned to the molecule and are used to determine the predicted toxicity value, which is obtained by averaging their toxicities. As explained in the discussion below, in both cases, the TPF finally selected might differ from that assigned during preassignment.
Results and Discussion The modeling strategy described was first applied to the training set to extract TPFs, build the priority matrix, and check its reliability. After the preassignment procedure (stage 1) on 240 compounds, Ntpf ) 133 different TPFs were identified; 123 fragments had one ACU, 8 had two ACUs, while 2 had 3 ACUs. For the sake of clarity, TPFs will be indicated here with the constituent ACUs in parentheses. For example, (209) refers to a TPF made up of the ACU 209, while (10, 3, 5) indicates a TPF made up of the ACUs 10, 3, and 5. Grouping many ACUs into one TPF had no important consequences as far as the training set molecules were concerned, since these ACUs had the same statistical distribution within the training set. It might, however, have significant consequences for the test molecules, as we will see later. To one training set molecule, missing a statistically significant ACU, no TPF was assigned; this compound (CAS No. 556-61-6) was therefore discarded from the training set, so the total number of training set molecules was 239. Once priorities were computed and retrieved in the priority matrix (accounting for 133 × 133 elements), we submitted training set molecules to postassignment (see stage 3 above). TPFs coming from postassignment were compared with those obtained during the preassignment (see stage 1). The comparison of pre- and postassignments enabled us to check the accuracy of the priority scheme from a qualitative point of view. In principle, the priority scheme should exactly reproduce the initial assignment. In practice, however, one might reasonably expect some molecules to show two different assignments, as the priority matrix does not account for individual molecular situations. Preassignment is a local procedure, in the sense that it deals with molecules as individual entities while assigning the TPF. In passing from the preassignment stage to priorities retrieval, we lose molecular identities because priorities are “averaged”
CAS No.
TPFPRE
TPFPOST
107534-96-3 1912-24-9 144-49-0 39515-41-8
(88) (7) (101) (50)
(110) (9), (207) (78), (101) (489, 490), (50)
on the whole training set. This operation is analogous to computing an integration over a finite domain. The result of such a mathematical operation depends on the limits over which the integral is computed. In the same way, the content of the priority matrix depends on the molecular domain, i.e., the training set. The priority matrix, therefore, has a global character instead of a local one, because it deals with the training set as a whole. The entire process is not reversible, because the local information about the role of each TPF in each molecule has been lost in building the priority matrix. For this reason, backtracking to the original assignment is not always possible, and postassignment and preassignment might be different. Applying our priority scheme to the training set, we found that postassignment matched preassignment for 186 out of 240 molecules. The accuracy, expressed as a percentage of correct reassignments, was therefore 77.5%; 54 molecules were erroneously reassigned. For the sake of discussion, 4 molecules not correctly assigned are reported in Table 2. In a comparison of the last two columns, where TPF reference integer numbers are reported, two different classes of erroneous reassignments come to light. The first contains molecules for which the assigned TPFs were completely different (for example, the first two). The second contains molecules for which pre- and postassignments partially agreed. The third and fourth molecules provide meaningful examples. In both molecules, two TPFs with the highest priority were found and assigned (see scenario “C” in postassignment), one correct and the other incorrect. The first type of incorrect reassignments is related to the presence of TPFs winning over many others. The high number of winning relationships rewards these fragments at the expense of all the others, finally spoiling the postassignment. Erroneous reassignments of the second type occur when the information retrieved in the priority matrix are not detailed enough to discriminate accurately TPFs. This especially happens when the proportion of conflicting or void relationships is higher than the winning or losing ones. The overall occurrence of these errors depends on the number and the structural diversity of the compounds under investigation. Qualitative assessment of the priority scheme has so far provided useful information about the causes of incorrect assignments. Nevertheless, we still have to check their effect on the predictive accuracy of the proposed QSAR. We want to quantify how incorrect reassignments affect the predicted toxicity values. We considered the models associated with preand postassignments. These models were built using average TPF toxicities as predicted molecular toxicities. We called the model preassignment (PRE) a “maximum performance model”, since the TPFs assigned were those having average toxicities closest to the experimental values. The PRE model therefore provided a solid reference for measuring postassignment
Top-Priority Fragment QSAR Approach
Chem. Res. Toxicol., Vol. 19, No. 11, 2006 1537
Table 3. Compounds with the Largest Absolute Differences between Observed (TOBS) and Predicted (TPOST) Toxicitiesa
a
CAS No.
TOBS
TPRE
TPOST
115-32-2 120-36-5 122-34-9 12683-31-7 14339-08-9 1563-66-2 1918-16-7 20018-09-1 2303-17-5 28249-77-6 58-36-6 82-66-6
3.475 2.672 0.456 2.353 3.217 2.765 3.095 3.511 2.405 2.390 5.157 2.082
3.543 2.718 1.394 2.367 3.187 2.718 2.833 2.760 2.382 2.382 4.524 2.382
2.137 0.867 1.605 3.387 4.601 4.053 1.869 2.273 4.449 1.215 3.688 3.403
Toxicity values coming from preassignment are also reported (TPRE).
performances. The correlation coefficient (R2) for this model was R2TRPRE ) 0.91. Interestingly, for the model obtained following the rules set for postassignment (POST), the result was comparable, R2TRPOST ) 0.85. This indicates that the overall effect of incorrect assignments on predicted values was small. The outcome was different when we restricted the comparison to the subset of the 54 incorrectly assigned molecules. In this case, we found a larger gap between pre- (R2SUBPRE ) 0.95) and post- (R2SUBPOST ) 0.63) assignment. A total of 12 compounds had an absolute difference between the observed and predicted (from postassignment) values larger than one LC50-96 h logarithmic unit. These are reported in Table 3. Most (8 of 12) of the incorrect postassignments reported in Table 3 belonged to the first class. For these compounds initial and final average TPF toxicities were very different. This suggested our approach was, in some cases, very sensitive to incorrect assignments of the first class, from a quantitative point of view. At the same time, since the number of compounds with large errors was very small (12) in comparison to that of training molecules (239), we concluded that the error due to first-class incorrect assignments was small. The results for preand postassignments can be found in Appendix A, where the name, CAS number, and observed toxicity have also been reported for each molecule. Given this result, we passed to the analysis of the test set. As specified above, 42 molecules were considered as the external test set and used for validating the model’s performance. These molecules were directly submitted to stage 3 of the analysis, where the search for TPFs was preliminarily carried out. For TPFs made up of more than one ACU, the presence of all constituent units was a necessary condition to assess their presence. For example, the TPF (10, 3, 5) is not present in molecules containing only the ACUs 10 and 3 or 3 and 5. Recognition of TPFs was successful for 41 molecules out of 42; the molecule for which no TPF was found (CAS No. 10706-2) could not be predicted. For five of the remaining compounds the assignment was multiple (see stage 3, scenario C) since two different TPFs were assigned. Only one TPF was associated with all other molecules. Among the TPFs assigned, one was made by two ACUs. According to stage 3, toxicity was then computed considering the TPFs with the highest priority. Predicted and observed toxicities for test molecules are reported in Appendix A for all compounds. Using these data and excluding the unpredicted molecule, the correlation coefficient was R2TSPOST ) 0.75. This result can be compared with that recently published by Toropov et al. (28). This study is, to our knowledge, the only one in the literature on this dataset and one of the few
Table 4. Results for Toropov et Al. and Our Model model
NTR
R2TR
NTS
R2TS
TPF Toropov et al.
239 233
0.85 0.77
41 41
0.75 0.64
applications of fragment-based method to pesticide risk assessment. The method proposed by Toropov et al. exploited the correlation weights of local and global graph invariants to predict toxicity. The size, in terms of number of molecules, of training and test sets was equal to that investigated in this work. Table 4 summarizes the results from both models for training and test sets in terms of correlation coefficients (R2TR and R2TS, respectively). The effective number of molecules predicted (NTR, and NTS) is also reported. Our model gave a better correlation coefficients for both training and test sets. This comparison, however, is more qualitative than quantitative. Despite both methods handle substructural features for modeling purposes, there is a huge difference between the methodological approaches adopted. While Toropov et al. approach simultaneously exploits all fragments for predicting the biological activity, our method makes use of one fragment to address the same issue. Looking at the results we obtained for this dataset, we concluded that an approach based on prioritization could be very effective in modeling pesticides, and our strategy is especially well-suited when toxicity is determined by the reactivity of a specific substructure and scarcely affected by the rest of the molecule. Even if both the configuration of the data set and our method do not strictly require a further statistical analysis, we performed a leave-one-out test. In this experiment all the compound in the training set were used to determine the TPFs and the priority scheme; the left-out compound was then predicted. The result is partly unsatisfactory (R2 ) 0.4) and, at a first sight, seems a sign of a scarce robustness of the approach. But, further looking at the worst predicted toxicities, we found that 23 out of the 29 compounds that lie far from the correlation line show a wrong prediction for the same reason: their removal from the training set causes the loss of the correct TPF from the analysis. This is again due to the high variability of the compound set, where some subclasses contain only two molecules. The consequence is that moving one of the two molecules to the test set causes the decrease of one ACU occurrences to one and, therefore, its elimination from the analysis (see above); the result is a wrong TPF selection. Five compounds are wrongly predicted in both the analyses, and only one compound is here badly predicted without a clear reason. If these 29 molecules are discarded, the LOO analysis shows an acceptable predictive potential (R2 ) 0.65). To properly judge our result, it is worth noting the factors that made modeling this dataset a real challenge. More than 20 pesticide classes are represented in the current dataset. These include organotins, organochlorines, organophosphates, carbamates, formamidines, terpenes, pyrethroids, phenols, spinosyns, pyrroles, pyridazinones, benzoylureas, etc. Some of those classes contain only one or two compounds. This clearly hampers the straightforward application of methods that demand structurally descriptive compounds to work properly, such as fragment-based ones. Considering these chemical features and the results, the prioritization strategy we presented was successful in modeling pesticide toxicity. Use of small ACUs enabled us to partially overcome the problem of dealing with a wide variety of structurally diverse compounds. Although simple, setting and retrieval of mutual priorities were accurate enough to discriminate TPF roles in different molecules. The
1538 Chem. Res. Toxicol., Vol. 19, No. 11, 2006 Table 5. Sketches of Five TPFs (Reference Indexes Bold Typed)a
a TPFs average toxicities (Tave) are reported in the second column. Training and test molecules to which the TPFs were assigned are also reported (the corresponding observed toxicities are shown in parentheses).
fact that the majority of the reassignments corresponds to the original ones (77.5%) indicates that our procedure can, within the limits of our main assumption, identify the expected answer. It also shows that the information retrieved in the priority matrix was appropriately exploited during the postassignment stage. It is important to test the validity of our main assumption. Some examples facilitate discussion. We chose the five TPFs sketched in Table 5. TPFs 60, 101, and 229 are made by a unique ACU; in contrast, TPFs (617, 618) and (608, 609, 610) provide two good examples of combined TPFs. All the TPFs are selected in the post phase. TPF 60 is present in some of the pyrethroid compounds and selects them very efficiently; TPF (617, 618) is also present in pyrethroid compounds and in the compounds selected by TPF 60. However, the compounds containing TPF (617, 618) are more toxic than the other pyrethroids, to which the TPF 60 was assigned. From the reactivity perspective, while TPF 60 is more representative and it could have been used to select all these compounds, TPF (617, 618) enhanced toxicological activity, as shown by the compounds reported in Table 5. It is therefore reasonable to conclude TPF (617, 618) to behave as a modulator rather than a toxicophore. TPF 101 is typical of aliphatic acids. It is also present with other ACUs in more molecules, but the toxicity of the selected compounds is reproduced more precisely by this single ACU. TPF 229 is a representative example of a small ACU characteristic of a subset of compounds. In principle, it can be used to select phosphates, phosphonates, thiophosphates, etc.; using the molecular toxicity, however, the choice can be limited to a subset that contains only thio derivatives. Finally, combined TPF (608, 609, 610)
Casalegno et al.
is very specific; in fact, the selected compounds are very similar and their toxicity values very close. We adopted the approach “one TPF-one toxicity”. If we assume that the correlation between toxicity and toxicophore is valid, we can also expect that the selected TPF properly matches the toxicophore of the current compound. However, from the results we cannot affirm that our model always selects the toxicophore in a compound. What we did is to group molecules using fragments AND toxicity. This means that the molecules that have the same TPF also share the toxicity. In some cases, this is coherent with the selection of the toxicophore, as demonstrated for the TPFs above; in other cases the connection between TPF and toxicophore is less clear. In these last cases the model can only select by toxicity; i.e., compounds that have similar toxicity but lack a fragment directly associated with a known toxicophore are nevertheless grouped using a TPF as the most probable common substructure. Then, the established TPFs are used to assign a toxicity to all compounds that, by the same mechanism, fall into the same group. The “one fragment-one toxicity” relationship works much better than our DFG approach (data not shown) (26). Our experience in toxicity prediction proves that the mechanisms that determine the toxicity are not parts of a single model; i.e. there are mechanisms based on the whole structure (e.g., narcotics) and mechanisms that specifically interact with a biological target (e.g., electrophilic reactants). In principle, identification of the real mechanism for each compound should lead to defined models and accurate predictions. However, this is not the case in toxicity prediction, where the mode of action (MOA) of chemicals is not always known. The results concerning the test set, in line with those recently published by Toropov et al. (28), look promising and encourage us to further improve the method.
Conclusions A new modeling strategy based on the prioritization of fragments contribution to toxicity was developed and applied to 282 pesticides for predicting their toxicity toward the rainbow trout (Oncorhynchus mykiss). The results for the training and test sets confirm the method’s reliability and efficacy from the quantitative point of view. Analysis of the TPFs assigned during postassignment enabled us to rationalize toxicity at the substructural level for some compounds. This, however, could not be done successfully carried out for all pesticides, because of the compounds’ structural diversity. We are currently working to improve the prioritization scheme, to fully exploit all available chemical information and deal better with heterogeneous datasets. Nevertheless, more data are needed to understand the toxicological properties of pesticides. Acknowledgment. We gratefully acknowledge the EU project DEMETRA (Contract No. QLK5-CT-2002-00691) for having partially supported and funded the current study. Supporting Information Available: Trout data set for 282 pesticides, divided into training and test sets. This material is available free of charge via the Internet at http://pubs.acs.org.
References (1) Food Quality Protection Act (FQPA), Public Law 104-170, Aug 3, 1996. (2) Council Directive concerning the placing of plant protection products on the market (91/414/EEC), OJ L 230, 19.8.1991, p 1. (3) Overview of the Ecological Risk Assessment Process in the Office of Pesticide Programs, Endangered and Threatened Species Effects
Top-Priority Fragment QSAR Approach
(4) (5) (6) (7)
(8)
(9)
(10) (11) (12) (13) (14) (15)
(16)
Determinations. Office of Prevention, Pesticides and Toxic Substances, Office of Pesticide Programs, Jan 23, 2004, Washington, DC. Guidance Document on Aquatic Ecotoxicology, under Council Directive 91/414/EEC, SANCO/3268/2001 rev. 4 (final), 17 Oct 2002. Guidance Document on Terrestrial Ecotoxicology under council Directive 91/414/EEC, SANCO/10329/2002 rev 2 (final), 17 Oct 2002. Guidelines for Ecological Risk Assessment (published on May 14, 1998: Fed. Regist. 63 (93), 26846-26924). Technical Guidance Document in support of Commission directive 93/67/EEC for new notified substances and Commission regulation (EC) No 1488/94 on risk assessment for existing substances. Brussels, Office for Official publications of the European Communities, 1996, ISBN 92-827-8013-9. Guidance for identifying pesticide chemicals that have a common mechanism of toxicity, for use in assessing the cumulative toxic effects of pesticides. Source: science policy issues and guidance documents (www.epa.gov/pesticides/trac/science). Assessment of the environmental properties and effects of pesticides transformation products. Sinclair, C. J., and Boxall, A. B. A., Final Report. Cranfield Centre for EcoChemistry, Cranfield University, Silsoe, Beds MK45 4DT, U.K. (Dec 2002). Vighi, M., Garlanda, M. M., and Calamari, D. (1991) QSARs for toxicity of organophosphorous pesticides to Daphnia and honeybees. Sci. Total EnViron. 109/110, 605-622. Nendza, M., Dittrich, B., Wenzel, A., and Klein, W. (1991) Predictive QSAR models for estimating ecotoxic hazard of plant protecting agents: target and non-target. Sci. Total EnViron. 109/110, 527-535. Nendza, M. (1991) Predictive QSAR models estimating ecotoxic hazard of phenylureas: aquatic toxicity. Chemosphere 23, 497-506. Devillers, J. (2004) Prediction of mammalian toxicity of organophosphorus pesticides from QSTR modeling. SAR QSAR EnViron. Res. 15, 501-510. Devillers, J., Pham-Delegue, M. H., Decourtye, A., Budzinski, H., Cluzeau, S., and Maurin, G. (2002) Structure-toxicity modeling of pesticides to honey bees. SAR QSAR EnViron. Res. 13, 641-648. Sparks, T. C., Crouse, G. D., and Durst, G. (2001) Natural products as insecticides: the biology, biochemistry and quantitative structureactivity relationships of spinosyns and spinosoids. Pest. Manag. Sci. 57, 896-905. Devillers, J. (2001) A general QSAR model for predicting the acute toxicity of pesticides to Lepomis macrochirus. SAR QSAR EnViron. Res. 11, 397-417.
Chem. Res. Toxicol., Vol. 19, No. 11, 2006 1539 (17) Devillers, J., and Flatin, J. (2000) A general QSAR model for predicting the acute toxicity of pesticides to Oncorhynchus mykiss. SAR QSAR EnViron. Res. 11, 25-43. (18) Perez, G. M., Gonzalez, D. H., Molina, R. R., Cabrera, M. A., and Ramos de Armas, R. (2003) TOPS-MODE based QSARs derived from heterogeneous series of compounds. Applications to the design of new herbicides. J. Chem. Inf. Comput. Sci. 43, 1192-1199. (19) Toropov, A. A., and Benfenati, E. (2006) QSAR models for Daphnia toxicity of pesticides based on combinations of topological parameters of molecular structures. Bioorg. Med. Chem. 14, 2779-2788. (20) Toropov, A. A., and Benfenati, E. (2006) QSAR models of quail dietary toxicity based on the graph of atomic orbitals. Bioorg. Med. Chem. Lett. 16, 1941-1943. (21) Debnath, B., Gayen, S., Basu, A., Ghosh, B., Srikanth, K., and Jha, T. (2004) Quantitative structure-activity relationship study using refractotopological state atom index on some neonicotinoid insecticides. Bioorg. Med. Chem. 12, 6137-6145. (22) Pearl, G. M., Livingstone-Carr, S., and Durham, S. K. (2001) Integration of computational analysis as a sentinel tool in toxicological assessments. Curr. Top. Med. Chem. 1, 247-255. (23) Mekenyan, O., Dimitrov, S., Schmieder, P., and Veith, G. (2003) In silico modelling of hazard endpoints: current problems and perspectives. SAR QSAR EnViron. Res. 14, 361-371. (24) Development of Environmental Modules for Evaluation of Toxicity of pesticide Residues in Agriculture (DEMETRA), project supported by the European Commission, Contract No. QLK5-CT-2002-00691. (25) Roncaglioni, A., Benfenati, E., Boriani, E., and Clook, M. (2004) A protocol to select high quality datasets of ecotoxicity values for pesticides. J. EnViron. Sci. Heal. B39, 641-652. (26) Casalegno, M., Benfenati, E., and Sello, G. (2005) An automated group contribution method in predicting aquatic toxicity: the diatomic fragment approach. Chem. Res. Toxicol. 18, 740-746. (27) Adamson, G. W., Lynch, M. F., and Town, W. G. (1971) Analysis of structural characteristics of chemical compounds in a large computerbased file, II. Atom-centred fragments. J, Chem. Soc. C 3702-3706. (28) Toropov, A. A., and Benfenati, E. (2006) Correlation weighting of valence shells in QSAR analysis of toxicity. Bioorg. Med. Chem. 14, 3923-3928.
TX0601814