In Silico Prediction and Experimental Evaluation of Furanoheliangolide Sesquiterpene Lactones as Potent Agents against Trypanosoma brucei rhodesiense

Antiprotozoal activity of STLs.The bioactivity data of 47 sesquiterpene lactones (compounds 1 to 40, 50, 51, 58 to 61, and 64) against T. brucei rhodesiense, T. cruzi, L. donovani, and P. falciparum as well as cytotoxicity against L6 rat skeletal myoblasts were published previously (7, 8, 15). Their activity data were included in the present QSAR analysis but are reported only as supplemental information (see Table SA1 in the supplemental material). For the present investigation, the in vitro activity against the mentioned target organisms/cells of 12 further STLs (compounds 41 to 45, 49, 52, 55, 57, 62, 63 and 65; for structures, see Fig. 1; for biological data, see Table 1) was tested. The pseudoguaianolides helenalin 1 and helenalin acetate 2 remained the most active compounds against T. b. rhodesiense in our series, with IC₅₀s of 0.05 and 0.06 μM (5, 7). Three further pseudoguaianolides included in the present series, namely, arnifolin, 2,3-dihydroaromaticin, and peruvin (compounds 41 to 43), displayed only a low level of activity. All three have a saturated ketone structure in the cyclopentane ring. Their low activity is thus in accordance with our previous structure-activity relationship studies, where the importance of a cyclopentenone structure for activity was already pointed out (7).

• Open in new tab
• Download powerpoint

• Open in new tab
• Download powerpoint

FIG 1

Structures of STLs used in the present QSAR analysis. Compounds 1 to 40 were included in our earlier QSAR study. They are shown here for a complete overview.

The highest antitrypanosomal activity within the newly tested series of STLs was displayed by the germacranolide lipiferolide (= 8β-hydroxyparthenolide acetate [compound 62]), a constituent of the tulip tree, Liriodendron tulipifera (Magnoliaceae), which showed an IC₅₀ of 0.23 μM against T. b. rhodesiense. It is noteworthy that this compound, recently also reported independently by others to exhibit activity against Plasmodium falciparum (21), was significantly more active against T. b. rhodesiense than four structurally related 9β-hydroxyparthenolide esters (compounds 58 to 61) recently isolated from Inula montbretiana (Asteraceae) (8) and also outmatched unsubstituted parthenolide (compound 34; IC₅₀ = 0.39 [7]).

QSAR analysis for activity against T. brucei rhodesiense.An attempt to test our previous quantitative structure-activity relationship (QSAR) model for compounds 1 to 40 (7) for external predictivity by calculating the activity data of 19 newly tested compounds (41 to 45, 49 to 52, 55, 57 to 65) failed to yield a significant correlation between experimental and predicted data (correlation coefficient [R] = −0.05; root mean square error of predictions [RMSEP] = 1.35 log units in pIC₅₀). In order to arrive at a more predictive QSAR model to be used for external predictions of a large set of compounds, an extended QSAR analysis was performed after the newly determined in vitro activity data were combined with those of our earlier studies (7, 8, 15) and data for 10 further STLs determined at the STPH laboratory under identical conditions in cooperation with other coauthors (compounds 46 to 48, 53, 54, and 56 [9, 10] and 66 to 69 [11]). Thereby, a much broader chemical diversity of the data set than that in our previous study was obtained. For complete data on activity against T. b. rhodesiense of all 69 STL, see Table SA1 in the supplemental material.

Thus, a coherent data set for 69 STLs was available for the present QSAR study. Three-dimensional molecular models for all compounds were generated and 123 molecular descriptors (for the full list, see Table SA2 in the supplemental material) were calculated from these 3D structures (see Materials and Methods). The set of descriptors used in our previous QSAR study (7) was extended with a series of 78 vsurf descriptors implemented in newer releases of MOE, which represent hydrophobic/hydrophilic surface properties calculated from molecular interaction fields as originally formulated by Cruciani et al. (22).

The data set was separated into a training set (n = 46 compounds) and a test set (n = 23). In order to ensure a maximum chemical diversity of the training and test set, this division was performed after diversity ranking based on molecular fingerprints as implemented in MOE (for the selected fingerprint scheme and other details, see Materials and Methods). The 46 most diverse compounds were chosen as the training set, leaving 23 molecules for a test set with similar chemical diversity. Subsequently, QSAR analysis was performed with the training set by applying a genetic algorithm/multiple linear regression (GA/MLR) approach (see Materials and Methods). For the present study, the GA was used to select a predefined number of descriptors (4, 5, or 6) from the whole matrix, that yield an optimal correlation (coefficient of determination [R²]) of the calculated and experimental data for the training set. The result of a GA/MLR calculation is a family or population, of predefined size, of QSAR equations ranked by their R² values. This population of optimized regression models is then submitted to leave-one-out cross validation (LOO-CV), yielding a cross validated coefficient of determination (Q² value) for each equation. In the present study, a significant increase in Q² was observed when the number of variables was increased from 4 to 5 (Q²_max rose from 0.567 to 0.658), whereas a further increase in the number of descriptors led to only a marginal increase, i.e., Q²_max = 0.667. Activity predictions for the test set molecules were then performed with the 10 best 5-variable equations. The best combination of internal and external predictivity was found for the following QSAR equation: pIC50(Tbr)=−3.630(±1.39)+0.649(±0.099)ENONCS+0.362(±0.091)ASAP4+0.643(±0.105)FASA−+0.670(±0.146)vsurf_D6−0.875(±0.159)vsurf_D8(1)

where Tbr is T. b. rhodesiense, n = 46, R² = 0.75, RMSE = 0.42, Q² = 0.65, and RMSEP = 0.50. Data were standardized to unit variance.

In this equation, the descriptor ENONCS refers to the accessible surface area of reactive carbon atoms in conjugated enone systems and thus represents a measure of the potential of an STL to react with nucleophilic groups of biological target molecules. ASAP4 is the fractional accessible surface area (23) attributable to atoms with a partial charge between +0.15 and +0.2 e (electrons). This charge interval is populated by hydrogen atoms attached to the double-bond carbons in α,β-unsaturated ketone groups such as H-2 and H-3 in helenalin (compound 1). This descriptor, already found to be influential in our previous antitrypanosomal QSAR (7), is also related to the presence of reactive structure elements. FASA⁻ is the ratio of accessible surface area contributed by atoms with negative partial charge in relation to the total accessible surface area. The positive regression coefficient indicates an enhancing effect of negatively charged surface on activity. Since in this series of molecules, negative partial charge is mostly caused by oxygen, this could be interpreted in terms of the propensity to accept H bonds. In any case, this is in agreement with our previous observation that descriptors of positive surface area had a detrimental effect (7). The descriptors vsurf_D6 and vsurf_D8 refer to properties of molecular interaction fields with a hydrophobic DRY probe (22). They correspond to the tendency to generate attractive hydrophobic interactions. In this context, the descriptor D6 corresponds to the volume of regions with an interaction energy of −1.2 kcal/mol with the DRY probe, whereas D8 is calculated at −1.6 kcal/mol. The two descriptors show some colinearity (R = 0.76), as can be expected from the fact that they reflect very similar molecular properties (for a matrix of pairwise correlation between all variables in the model, see Table SA3 in the supplemental material). However, it was found that neither of them alone led to a significant improvement in comparison with a 3-variable model consisting only of the descriptors already discussed but that both have to be present in the equation to yield a very significant increase of the coefficients of determination. Their t test P values (both ≪10⁻⁴ [see Table SA4 in the supplemental material]) show that they both make highly significant contributions to the QSAR equation's overall performance. Although this combination of two descriptors is somewhat difficult to interpret, the positive and negative signs of their regression coefficients and the fact that they are both necessary for the model might indicate that it is more favorable for activity if a molecule has a relatively large region causing hydrophobic interactions of medium strength (vsurf_D6) than a large area forming strong hydrophobic interactions (vsurf_D8), i.e., that the lipophilicity should not be too high.

A plot of the experimental versus calculated pIC₅₀ values for the training and test sets is shown in Fig. 2 (the underlying data are reported in Table SA2 in the supplemental material). As can be seen, the model performs satisfactorily to predict the activity of most test set molecules. The predictions are quite reasonable, especially in the region of molecules with high activity (i.e., a pIC₅₀ of >6). Overall, however, the coefficient of determination, P², for the predictions was only 0.35, which was due to three out of 23 compounds (no. 7, 32, and 56) whose predictions were of rather poor quality (i.e., more than 1 log unit too low). The model's failure to predict the activity of compounds 32 and 56 correctly can readily be attributed to the fact that their structural properties are not well represented in the training set. Thus, the predicted activity of compound 32 is much too low, since this compound differs only in stereochemistry from the much less active compound 31 (1-α-OH versus 1-β-OH), a difference which is not accounted for by the descriptors of the QSAR model. Compound 56 contains an ester side chain (dihydroxytigloyl) which does not occur in any other compound of the set but may be responsible for its increased activity in comparison with some similar compounds. In case of compound 7, the poor prediction is not as straightforwardly explained, since its congeners in the training set (compounds 5 and 6 with respect to the STL core and compound 4 with respect to the ester side chain) represent its structure quite well. Their calculated activity values are, however, generally somewhat too low, so that several factors leading to an underestimation of activity may sum up in this case.

• Open in new tab
• Download powerpoint

FIG 2

Correlation plot of experimental versus predicted pIC₅₀ values for anti-T. brucei rhodesiense activity of 69 STLs as calculated by QSAR equation 1. Diamonds and circles, model calibration and leave-one-out (LOO) predictions, respectively, for the training set (46 compounds); triangles, predictions for the external test set (23 compounds); squares, three test set compounds (compounds 7, 32, and 56) with particularly poor predictions. Omission of these outliers yielded a much better correlation coefficient for the remaining 20 compounds [test set 2].

After omission of these outliers, P² rose to 0.62 and was thus compatible with the Q² value of the internal cross validation predictions.

It therefore appeared reasonable to attempt predictions for yet-untested STLs in order to find new candidates with high predicted activity for testing.

Activity predictions for a virtual library of 1,750 STL structures.Our in-house database of STLs, originally created for the purpose of structure elucidation, contains optimized 3D molecular structures for a collection of almost 2,100 STLs published by Budesinsky and Saman along with their ¹³C NMR data (19). For the present study, structures violating Lipinski's rule of 5 (20) were eliminated from the database, so that 1,750 structures remained for further consideration. Using the QSAR model defined by equation 1 above, activity predictions were made for these “druglike” STL structures.

On the whole, 71 compounds were predicted to have a pIC₅₀ value against T. b. rhodesiense equal to or larger than 7.0 (i.e., IC₅₀s of ≤0.1 μM). Quite notably, among these 71 structures, 15 were found to belong to the furanoheliangolide subgroup of germacranolide-type STLs, of which, to the best of our knowledge, no representatives had hitherto been tested for activity against T. b. rhodesiense.

On this background, four furanoheliangolides (compounds 70 to 73 [Fig. 1]) available in sufficient purity and quantity were selected for in vitro testing against this parasite.

Experimental confirmation of activity predictions for furanoheliangolide-type STLs.The furanoheliangolide goyazensolide (compound 70), with a predicted pIC₅₀ of 7.62 (i.e., IC₅₀^pred = 0.02 μM) was chosen for testing along with its structural congener budlein A (compound 71), which was predicted to be less potent (pIC₅₀^pred = 6.60; IC₅₀^pred = 0.25 μM) but still active enough to make it interesting to test. Additionally, a structural analogue of compound 70 in which the potentially reactive α,β,γ,δ-unsaturated ketone system is reduced to an α,β-unsaturated system, i.e., 4,5-dihydro-2′,3′-epoxy-15-deoxygoyazensolide (compound 72), was also chosen. Its activity was predicted to be considerably lower (pIC₅₀^pred = 6.17, i.e., IC₅₀^pred = 0.68 μM) due to the hydrogenation of the conjugated γ,δ-double bond between C-4 and C-5. Furthermore, a structural analog with an exocyclic orientation of the γ,δ double bond (i.e., Δ^4,15 instead of Δ^4,5 as in compounds 70 and 71), namely, 4,15-isoatripicolide tiglate 73 (pIC₅₀^pred = 7.04; i.e., IC₅₀^pred = 0.09 μM), was also tested.

In the T. b. rhodesiense bioassay, compounds 70, 71, and 73 were highly active, with IC₅₀s of 0.073 and 0.072 μM, respectively (70 and 71; pIC₅₀ = 7.14) and of 0.015 μM (73; pIC₅₀ = 7.82). Comparison with the predicted values shows that the QSAR model slightly overestimated the activity of compound 70 and underestimated that of compounds 71 and 73. However, all three had correctly been predicted to show very high activities comparable to those of compounds in the upper range of the training set, so that the validity of the QSAR model was also confirmed experimentally with these truly external predictions. At the same time, the 4,5-dihydrogoyazensolide derivative compound 72 was considerably less active than compounds 70 and 71, as predicted. Although its measured IC₅₀ of 0.20 μM was somewhat lower than the predicted one, this prediction once more confirms that the chosen QSAR model allows external predictions with reasonable accuracy.

Comparison of compound 70 with compound 71 shows that the nature of the acyl moiety as well as the stereochemistry at C-8 does not have any significant influence on anti-T. b. rhodesiense activity, whereas changes in the α,β,γ,δ-unsaturated carbonyl structure as observed in compounds 72 and 73 had a higher impact. In case of compound 72, the lower activity is readily explained by the loss of one potentially reactive site at C-5. Better accessibility of the reactive exocyclic methylene group C-15 for Michael addition to a potential trypanosomal target can be held responsible for the increased activity of compound 73 in comparison with the activated methine carbon (C-5) in compounds 70 and 71.

Compound 73 is the STL with the highest in vitro activity against T. b. rhodesiense identified until now. It also displays a much higher selectivity index (SI = 77) than its congeners 70 and 71 (SI = 7 and 10, respectively) as well as the most active pseudoguaianolides reported so far (e.g., helenalin 1 and helenalin acetate 2, with SIs of 19 and 13, respectively).