Computational Prediction of Drug Solubility in Fasted Simulated and Aspirated Human Intestinal Fluid

Purpose To develop predictive models of apparent solubility (Sapp) of lipophilic drugs in fasted state simulated intestinal fluid (FaSSIF) and aspirated human intestinal fluid (HIF). Methods Measured Sapp values in FaSSIF, HIF and phosphate buffer pH 6.5 (PhBpH6.5) for 86 lipophilic drugs were compiled and divided into training (Tr) and test (Te) sets. Projection to latent structure (PLS) models were developed through variable selection of calculated molecular descriptors. Experimentally determined properties were included to investigate their contribution to the predictions. Results Modest relationships between Sapp in PhBpH6.5 and FaSSIF (R2 = 0.61) or HIF (R2 = 0.62) were found. As expected, there was a stronger correlation obtained between FaSSIF and HIF (R2 = 0.78). Computational models were developed using calculated descriptors alone (FaSSIF, R2 = 0.69 and RMSEte of 0.77; HIF, R2 = 0.84 and RMSEte of 0.81). Accuracy improved when solubility in PhBpH6.5 was added as a descriptor (FaSSIF, R2 = 0.76 and RMSETe of 0.65; HIF, R2 = 0.86 and RMSETe of 0.69), whereas no improvement was seen when melting point (Tm) or logDpH 6.5 were included in the models. Conclusion Computational models were developed, that reliably predicted Sapp of lipophilic compounds in intestinal fluid, from molecular structures alone. If experimentally determined pH-dependent solubility values were available, this further improved the accuracy of the predictions.


INTRODUCTION
Modern drug discovery programs using high throughput screening and combinatorial chemistry continue to favor the selection of large and lipophilic new chemical entities (NCE). This is in spite of their poor aqueous solubility [1][2][3][4], the increased awareness of related problems and the multitude of mnemonic rules for avoiding these compounds with low or variable absorption and pharmacokinetics [5][6][7]. Solubility in intestinal fluids is a key property for estimating absorption of oral drugs and in this context aspirated human intestinal fluid (HIF) has been described as the gold standard medium for these estimations [8,9]. However, there are disadvantages associated with its use. Ethical concerns about HIF sampling regulate the availability of the fluid. The amount available to the scientific community is therefore sparse and expensive. Other concerns are its low buffer capacity and batch variation in pH and bile content. These variations are due to differences in aspiration protocols and individual differences between volunteers [9,10], although differences due to the latter can be somewhat alleviated by pooling samples. Alternatives to HIF for dissolution testing and solubility measurements became available when Dressman and colleagues introduced fasted state simulated intestinal fluid (FaSSIF) in 1998 [11]. Other biorelevant dissolution media (BDM) have been developed since then, some to closer mimic the intestinal milieu [12,13] and others for ease of preparation or lower expense [14,15]. Nonetheless the use of FaSSIF continues to be widespread and a large number of compounds and formulations have been evaluated in them. The medium contains taurocholate and lecithin that form mixed lipid aggregates in the form of vesicles which are colloidal structures known to efficiently solubilize drug molecules. The extent of solubilization is dependent on lipid concentration and substance-specific properties such as size, charge, flexibility, and lipophilicity. Drugs with a partition coefficient between octanol and water (logP) greater than 3 have considerably higher apparent solubility (S app ; the total concentration of drug dissolved in the lipid-containing dissolution medium) in BDM than in water or buffers [16][17][18][19]. Solubilization in lipid aggregates and molecular ionization as a response to the pH of the fluid can increase the S app in BDM of compounds several orders of magnitude compared to that observed in water.
Solubility measurements are time consuming. More importantly, the substances of interest must be synthesized before their solubility can be evaluated. Computational predictions of solubility, on the other hand, are rapid and can be performed on large compound libraries without synthesis of the substances. This provides the medicinal chemists with solubility profiles on which they can make better informed decisions, and the costs associated with the pharmaceutical profiling cycle are reduced because of decreased demand for expensive simulated or aspirated intestinal fluid. Solubility predictions in BDM, or even better HIF, are therefore highly warranted.
Numerous models for the prediction of intrinsic aqueous solubility (S 0 ), i.e., the solubility of the neutral compound, have been developed [20]. One of the most renowned is the general solubility equation (GSE) [21] and derivatives thereof [22,23]. These are based on logP, which typically can be computationally predicted with an RMSE of around one log 10 unit, and experimentally determined melting temperature (Tm). The latter is sometimes replaced with other properties more amenable to trivial calculations or prediction such as MW [24]. Aqueous solubility can also be successfully predicted using calculated molecular descriptors, see e.g. [25,26],. However, in the gastrointestinal tract, pH values range from 2.5 in the stomach to~6.9 in the jejunum [10]. This pHgradient greatly impacts the ionization of protolytic compounds and hence, the observed S app is dependent on the extent of ionization of a particular molecule. The pHdependent solubility can be calculated from S 0 and the dissociation constant (pKa) with the Henderson-Hasselbalch equation [27]. However, the accuracy of these estimations varies considerably because the Henderson-Hasselbalch equation does not take into account aggregation or common ion/salt effects [28]. The complexity increases even more when solubility is measured in BDM since the apparent solubility in these media is a result of ionization, aggregation and solubilization. We have previously attempted to predict S app in biorelevant media using a small dataset [16,17]. A predictive artificial neural network (ANN) model for FaSSIF solubility is also available in the commercial software ADMET Predictor from Simulations Plus. However, no transparent models for prediction of FaSSIF S app have been developed using publicly available solubility data for drugs, nor have any predictive models of HIF S app been proposed. Here, we report an open database applicable for solubility modeling in FaSSIF and HIF. This database has been used to develop transparent and reliable models for the prediction of solubility in FaSSIF and HIF with the aim of revealing molecular features that drive solubilization in these fluids.

Datasets
S app values for 86 drugs in FaSSIF (3 mM taurocholate, 0.75 mM lecithin in PhB pH6.5 [11]) were extracted from inhouse databases [8,16,17,29] and literature sources [30][31][32][33][34][35][36][37][38] ( Table I). To reduce experimental variability in the dataset the main part of the compounds was obtained from our in-house databases in which solubility measurements taking use of shake-flask or the μDISS Profiler are reported. Only compounds with a calculated logP greater than 2 were included since it is assumed that there is significant solubilization of highly lipophilic compounds in the mixed lipid aggregates present in FaSSIF [16][17][18][19]. Hence, we argue that for compounds with log P<2, in silico models predicting solubility in pH-adjusted water/simple buffer are also predictive of their solubility in intestinal fluid (Fig. S1). The FaSSIF S app values were supplemented with the corresponding S app values in PhB pH6.5 and the Tm, when available, for the free base or free acid (i.e., not salts).
This same approach was used in the selection of a dataset for which solubility in aspirated HIF was available (Table I). An additional criterion for this dataset was that corresponding solubility measurements were available in FaSSIF and that the same protocol had been used for solubility measurements in FaSSIF and HIF. The final HIF dataset consisted of 48 compounds. In addition to the FaSSIF and HIF literature datasets, a discovery dataset of 26 AstraZeneca proprietary compounds for which the S app in FaSSIF, HIF and PhB pH6.5 were available was used as an external validation dataset for the developed models.

Calculation of Molecular Descriptors
Molecular structures for the compounds were obtained as SMILES strings (Table SI) and converted to energy-minimized, three-dimensional structures with added implicit hydrogens using Corina 3.49 (Molecular Networks, Germany). Molecular descriptors from the resulting structures were generated using DragonX 6.0.16 (Talete, Italy). The descriptors were blinded to avoid selection bias, cubic-root transformed, mean-centered and scaled to unit variance. Thereafter, any descriptors not displaying normal distribution were excluded from the model development.    We also used ADMET Predictor 6.5 (Simulations Plus, CA) to predict FaSSIF S app for the compounds and to calculate pKa (Table SI) and pH-dependent lipophilicity at pH 6.5 (logD pH6.5 ), polar surface area (PSA), and number of rotatable bonds, see Table I. The predictions of FaSSIF were used as a comparator for our predictions.

In Silico Model Development
Partial least squares projection to latent structures (PLS) models were developed with the purpose of predicting S app values in FaSSIF and in HIF. The respective models were developed with Simca-P 13.0.2.0 (Umetrics, Sweden) using a standardized protocol previously implemented by our group [39,40]. The responses were used in the logarithmic form of the solubility in the two different media. For the FaSSIF model, the compounds were randomized into a training (Tr, n=56) and a test (Te, n=30) set. The compounds were sorted by their S app to achieve an even distribution and as wide a predictive range as possible. Every third compound was then assigned to the test set. The structural diversity and the suitability of the selected test set was thereafter tested by a principal component analysis (PCA) extracted from all descriptors. Any training set outliers identified in the PCA (Supporting Information, Fig. S2) or distance-to-the-model-of-X [41] (DmodX) were moved to the corresponding test set to avoid such compounds weighting the training of the models. PCA was also used to ensure that the training and test sets were well distributed in the chemical space. In addition to the literaturederived test set, a discovery dataset of 26 AstraZeneca inhouse compounds was used to challenge the model. For the HIF model, these discovery compounds were used as the sole test set while the training set consisted of the literature values. PCA confirmed that the chemical space of this dataset was covered by the training set (Fig. S2).
All but the top 100 descriptors (Table SII) (as defined by variable of importance-to-projection, VIP) were excluded in the first step of model development. The variable selection procedure was thereafter based on the VIP and the loading plots, and monitored by the leave-one-out (using 7 groups), cross-validated R 2 (Q 2 ). The variable selection was performed to remove non-significant descriptors and increase model robustness. If the exclusion of a variable resulted in an equal or improved Q 2 , the variable was permanently eliminated from the model. The variable selection procedure was repeated until no further descriptors could be removed without a resultant lower Q 2 . Only thereafter was the accuracy of the prediction of the test set investigated.
In the second step, the impact of experimental data commonly available during early development was investigated by adding such measured data to the final model obtained after the completion of variable selection. The investigated experimental data were S app in PhB pH6.5 , Tm, and logD pH6.5 . Measured PhB pH6.5 S app was available for 76% of the compounds (Table I). It was not possible to extract measured logD pH6.5 for a large number of the compounds and therefore the calculated logD p H 6 . 5 (ADMET Predictor, Simulations Plus, CA) was used for all of them. Whether these properties were beneficial or not for the models was evaluated in the same manner as for the calculated descriptors, see above.
In addition to the models above, a consensus model was established which used the developed FaSSIF PLS model based on calculated descriptors only and the predictions obtained from ADMET Predictor ANN model. No weighting was performed. The consensus model used the average of the predicted logS app from the developed PLS model and the commercial ANN model.

Physicochemical Properties and Apparent Solubility
The datasets used were structurally diverse. The FaSSIF modelling dataset (Table I) ranged in size from 206.3 to 1202.6 Da with a median mass of 369; hydrogen bond capacity in the form of PSA ranged from 2.2 to 325 Å 2 with a median of 122 Å 2 ; and molecular flexibility (described by the number of rotatable bonds count) ranged from 0 to 22, with a median of 5. The dataset was selected to focus on lipophilic compounds and therefore all compounds had a calculated AlogP-value >2 (obtained from the software DragonX). In spite of this lipophilicity criterion, the predicted pH-dependent  (Table I).
The solubility range of the training set was similar to that of the literature test set, whereas the 26 discovery compounds used to challenge the models had a somewhat lower solubility (Fig. 1).
The solubility values of the dataset revealed modest relationships between measured solubility in PhB pH6.5 and FaSSIF or HIF (R 2 :0.61 and R 2 :0.62, respectively, Fig. 2). The PhB pH6.5 and FaSSIF S app correlation was considerably weaker for the more lipophilic compounds (logD pH6.5 >4) for which R 2 decreased to 0.28. However, for compounds with a logD pH6.5 <3, there was a strong correlation (R 2 : 0.82) The correlation between FaSSIF S app and HIF S app was also strong (R 2 : 0.78) but unaffected by lipophilicity (Fig. 2c). Under-or over-prediction of the solubility was not related to melting point or lipophilicity.

Prediction of Intestinal Solubility
The developed PLS models are summarized in Table II. The FaSSIF model (Fig. 3a) required seven calculated descriptors to produce two principal components resulting in R 2 of 0.69, Q 2 of 0.64 and an RMSE tr of 0.48 log 10 units. The HIF model (Fig. 4a), based on nine descriptors, had a higher R 2 of 0.84 (Q 2 of 0.78) and a lower RMSE tr (0.34). The inclusion of experimentally determined S app in PhB pH6.5 strengthened the predictive power of both models to R 2 of 0.76 (FaSSIF) and 0.86 for HIF (Figs. 3b and 4b, respectively), and it reduced the RMSE of the test sets (Table II). Inclusion of Tm and logD pH6.5 did not improve the developed models. Further, these properties were unable to identify over-or under-predicted compounds or any clusters.
All descriptors remaining after the variable selection were significant in either both or the last component. For the FaSSIF model these include: i) eigenvalues weighted by bond order (Eig04_EA(bo)) or edge degree (Eig04_AEA(ed)); ii) a spectral moment of order 6 from Burden matrix weighted by van der Waals volume; and iii) a second-component accessibility directional WHIM index weighted by van der Waals volume; all these descriptors negatively influenced the solubility. These descriptors are related to some extent to molecular size. Geary autocorrelation (GATS4s) and Morse signal  Fig. 5a and b). The HIF model descriptors included an eigenvalue weighted by bond order (Eig03_AEA(bo)), Morse signal 4 weighted by ionization potential (Mor04i) and an eigenvalue from reciprocal squared distance matrix (SpMax_H2); these were found to limit the S app . The HIF S app was further positively influenced by the following descriptors:, R maximal index (RTe+) and a WHIM index (E3e) both weighted by Sanderson electronegativity, R maximal autocorrelation weighted by ionization potential (R1i+), the presence of N -O at topological distance 5(B05[N-O]), CATS2D acceptoracceptor at lag 09 (CATS2D_09_AA) and hydrogen attached to alpha carbons (H-051) (Fig. 5c and d).
The model using calculated descriptors over-predicted the FaSSIF S app of albendazole and cilostazole by one logarithmic unit or more. On the other hand, ivermectin, tamsulosin, and tolectin were all under-predicted by one log 10 unit or more. The inclusion of measured S app in PhB pH6.5 as a descriptor improved the predictions for cilostazole and albendazole and reduced their residual values by 0.41 and 0.69 log 10 units respectively.
The consensus model based on the calculations obtained from our FaSSIF model and the ADMET Predictor results exhibited a lower RMSE te (0.70) compared to each of the models separately (Table III). The PLS model was however more accurate compared to the consensus model in the prediction of the solubility of the discovery test set. Importantly, predictions of outliers resulting from each of the PLS and ANN models were greatly improved by consensus modelling (Fig. 6).

DISCUSSION
Intestinal solubility, together with permeability over the intestinal wall, are the two most important drug properties determining absorption after oral intake. Solubility measurements in HIF will continue to be important in understanding intestinal solubility, but the medium is expensive and subject to batch variations. Therefore, a number of BDMs have been developed as robust and reproducible surrogates. The strong correlation between FaSSIF and HIF S app found herein confirms those reported previously [8,9] and further supports the use of in vitro experiments in BDMs, such as FaSSIF, for the prediction of intestinal solubility.
The aim of this study was to develop predictive models for S app in FaSSIF and HIF using calculated descriptors alone or in conjunction with experimental data likely to be available in early drug discovery or development stages. The descriptors included in the final FaSSIF model can be used to interpret molecular properties of importance for solubility in FaSSIF. The descriptors reveal that larger structures are solubilized to a lesser extent than the smaller ones. Most likely this is as a result of the increased cavity that needs to be formed in the  water as well as the increased molecular surface area exposed to the water. In addition, aromatic structures were revealed to be less hydrated than aliphatic ones. We speculate that this could be because of their stronger crystal lattices, due to the stronger van der Waals interactions formed by the dense packing. Further, when in the water, rigid aromatic structures have a larger molecular surface area exposed to the solvent than flexible aliphatic chains that can change conformation to shield the carbon skeleton from water molecules. The descriptors also identify the importance of hydrogen bond donor and acceptors for the hydration of the molecule. Although the descriptors in the HIF model differ to some extent from those of the FaSSIF model, they too reflect similar properties.
During the model development we tested three empirical or semi-empirical descriptors. The addition of experimentally determined S app in PhB pH6.5 as a descriptor improved both developed PLS models considerably. This solvent can be regarded as a blank FaSSIF because it is a phosphate buffer (pH 6.5) that does not contain any mixed lipid aggregates. Ionizable compounds are therefore charged to the same extent in both media and information on this effect facilitates the prediction. The relationship between PhB pH6.5 and FaSSIF S app was modest and weak for the lipophilic compounds. It is therefore interesting to note that such a divergence was not seen in the PLS predictions regardless of whether or not PhB pH6.5 was included as an experimentally determined descriptor. Tm and lipophilicity expressed as logP or logD pH6.5 are common inputs for aqueous solubility predictions, and were therefore tested herein. However, neither logD pH6.5 nor Tm improved the developed models. Intestinal fluids contain lipid aggregates that may solubilize lipophilic drug molecules and it is well-known that logP is not a good descriptor of solubility in lipids, see e.g. [40],. The reduced solubility of the bulk water seen in intestinal fluids with higher logD pH6.5 is to some extent compensated by partitioning to and/or solubilization in the aggregates. To elucidate if Tm or logD pH6.5 were described to some extent by the selected calculated descriptors these properties were used as responses in the models. We found no indication of correlation between the selected descriptors and lipophilicity or solid state properties of the compounds (R 2 <0.45). An interesting aspect of the increased accuracy of the predictions when S app in PhB 6.5 is included as descriptor is that the influence of the solid state on the S app is diminished. Hence, the PhB 6.5 contributes to better predictions at two levels; the hydration is better described as a result of the correct description of the pH-dependent solubility and the impact of the dissociation of molecules from the crystal lattice is embedded in this solubility input variable.
It was possible to further improve the FaSSIF predictions by performing consensus modeling based on the combination of two different computational models that used only calculated descriptors. The combination of the (ANN) FaSSIF model results (obtained from ADMET Predictor) with the PLS predictions increased the predictive power, as identified from the lowered RMSE te (Table III). Three test sets were evaluated: all test compounds (n=75), a literature-derived test set (n=49), and a discovery test set of AstraZeneca proprietary compounds (n=26). For all three, the consensus model performed better than the worst-performing computational model. Most importantly, the consensus model increased the robustness and reduced the number of outliers (Fig. 6). Indeed, for the 23 compounds that were 10 to 320-fold over-or underpredicted by either of the two models, the employment of the consensus model reduced the RMSE to <1 log unit for 12 of them. Of these 23 compounds, 14 were significantly falsely predicted by the ANN, 5 by the PLS, and 4 by both. Since it is difficult to deem beforehand which one of several models will be the most accurate predictor for any new compound or compound series, it is advisable to employ consensus modelling based on two or more models.
There are a number of hurdles to allow increased accuracy in solubility predictions in BDM. The models developed herein that are based on calculated descriptors alone are reliable and statistically strong, but the observed residual errors do imply that the predictions could be up to tenfold off in either direction. This is not uncommon for predictions of solubility in aqueous media and the developed models are certainly accurate enough for guiding decision-making in drug discovery and development. To further improve the computational predictions of solubility in media containing mixed lipid aggregates molecular dynamics (MD) simulation is a promising tool to study solubilization interactions [42,43]. MD simulation also has the potential to predict selfaggregation and the tendency to form mixed micelle aggregates. Food effects on bioavailability is another issue related to BDM and larger datasets need to be studied for solubility in fed state BDM to obtain information on molecular features of drug molecules that are significantly affected by the increased lipid content in the fed state.

CONCLUSIONS
Measurements of solubility in physiologically relevant media such as FaSSIF and HIF are costly and in part restricted by the limited access to intestinal fluids aspirated from donors.
Herein we present computational approaches to instead predict intestinal solubility taking use of calculated molecular descriptors to allow this property to be estimated already before compound synthesis. It was found that the S app of lipophilic compounds in FaSSIF and in HIF was possible to predict by this approach. The most accurate predictions were obtained when a consensus modeling approach was used, which reduced the number of outliers obtained from predictions based on a single computational model. Further, we have examined experimental parameters within reach during early drug development and identified pH-dependent solubility as a descriptor that further increases the accuracy of the predictions. Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.