Chem-Bio Informatics Journal
Online ISSN : 1347-0442
Print ISSN : 1347-6297
ISSN-L : 1347-0442
Model Repository
Prediction Models for Fraction of Absorption and Membrane Permeability using Mordred Descriptors
Tsuyoshi EsakiTomoki YonezawaDaisuke YamazakiKazuyoshi Ikeda
著者情報
ジャーナル フリー HTML
電子付録

2022 年 22 巻 p. 46-54

詳細
Abstract

Machine learning (ML) models are cost-effective methods that have accelerated the identification of novel drug candidates in pharmaceutical research. These in silico methods estimate the characteristics of chemical compounds using calculated physicochemical features or using molecular sub-structure fingerprints. This rise in the deployment of machine learning models is facilitated by the development of numerous ML packages that enables researchers to build local models to meet their requirements. Despite the growing ease of building ML models, programming these informatics-driven solutions can be arduous for wet-lab researchers. In this study, we present a template for ML model construction that would enable researchers to efficiently reproduce ML models. We constructed prototype models to estimate the fraction of absorption and membrane permeability of a chemical compound using Mordred descriptors.

1. Introduction

Identifying potential drug candidates using cost-effective methods is the need of the hour in pharmaceutical research. In this regard, machine learning (ML) methods, including deep learning (DL), are an effective in silico drug discovery strategies to characterize chemical compounds using calculated features. The pharmacokinetic parameters related to human intestinal absorption, such as the fraction of absorption (Fa) and membrane permeability measured using the Caco-2 cell line (Papp), directly affect bioavailability and dosage schedules. Owing to this, early-stage estimation of these parameters using chemical features alone is desired.

Although ML models to predict pharmacokinetic parameters, such as Chemi-Net [1] and SwissADME [2], have been developed, they require the user to upload compound information to a web server. This pre-requisite invites concerns relating to confidentiality. Additionally, numerous packages to perform ML have supported an environment that fosters the construction of prediction models in each company. DL models have been created by winners of Kaggle contests [3] that require large measurement datasets [4]. Collecting such extensive experimental data in pharmaceutical research requires considerable effort and resources. Therefore, ML is given a higher priority for constructing models because of its ease of parameter tuning.

Additionally, these ML packages enable researchers to construct local in silico models to meet specific needs with the predictive capability for specific chemical compounds. These local models enable researchers to find promising candidates to target diseases for each company. However, deployment of these models requires expertise in programming in addition to performing trial-and-error for data preparation, descriptor selection, and hyperparameter tuning. These hurdles can prevent a wet-lab scientist from realizing the predictive potential of ML models.

To bridge this gap, we present a template for model construction that would enable any researcher to reproduce an ML model and interpret the results quickly [5]. Jupyter notebook [6] and descriptors have been shared as supporting information to benefit the researcher.

2. Materials and Methods

The data management and preprocessing were performed using the Ubuntu Linux operating system version 18.04 with Intel Xeon (3.9 GHz to 4.2 GHz) processor and 32 GB of RAM. We used previously published datasets [7] to construct two models to predict Fa and Papp of drug candidates. To perform molecular-descriptor calculation and model construction, we used Python (ver. 3.9) [8] in Jupyter notebooks. We used a commercial software to generate 3D conformer structures, owing to which the compound structure was not shared. However, the generated descriptors and scripts are made publicly available as supporting information. The details of model construction are briefly described in the following sections.

2.1. Overview of framework

In this study, two regression models were constructed (this has been shared in the supporting information file). Structure data format (SDF) is required as input file to reproduce our models. The calculated descriptors were exported in comma-separated values (CSV) format. Two regression models (Fa and Papp) were constructed using this CSV file, and proper descriptors and machine learning algorithms were decided.

Figure 1. Overflow of framework

The top two yellow nodes show the data collection and descriptor calculation. The left blue nodes and arrows show the workflow to construct Fa regression models. The right green nodes and arrows show the workflow to construct Papp regression models. The numbers of figures and tables in this overflow are linked in this article.

2.1. Libraries

In this study, six libraries were employed; NumPy (ver. 1.20.3), pandas (ver. 1.4.1), scikit-learn (ver. 1.0.2) [9], RDKit (ver. 2021.09.5) [10], Mordred (ver. 1.2.0) [11], and SHAP (ver. 0.40.0) [12].

2.2. Datasets

The dataset contained information on the chemical structures of 5567 compounds as SMILES strings. In this dataset, the number of Fa and Papp experimental values were 946 and 4460, respectively. While searching for tools for 3D conformer generation, CORINA (ver. 4.4.0) successfully generated the highest number of 3D descriptors using Mordred when compared with RDKit [13]. Owing to this result, we used CORINA to generate 3D structures of the chemical compounds as SDF.

2.3. Descriptor calculation

We used the molecular descriptor calculator, Mordred [11], to calculate descriptors (1D and 2D descriptors: 1613; 1D, 2D, and 3D descriptors: 1826). These descriptors were calculated for both Fa and Papp datasets.

2.4. Data preparation

The Fa dataset was randomly split into training (70%, 660 compounds) and test set (30%, 286 compounds) using the train_test_split function in scikit-learn. The Fa measurements ranged between 0.0 and 1.0, with localization around either 0.0 or 1.0. We transformed these values to log10(Fa/(1 - Fa)) to scatter the response variable. Additionally, Fa = 0.0 was set to 0.01, and Fa = 1.0 was set to 0.99. The distributions of this dataset in terms of molecular weight and SlogP that were calculated using Mordred are shown in Figure 1. The p-values of the Wilcoxon signed rank test were higher than 0.05 in the converted Fa, MW, and SlogP measurements, confirming that there was no bias in the training and test sets.

Figure 2. Distribution of Fa datasetThe x-axes in the figures above (left to right) show the Fa values, converted Fa, molecular weight, and SlogP

The y-axes indicate the number of samples. The blue and orange bars represent the training and test sets distribution, respectively.

The range of Papp measurements was between 0.00016 and 880 [10-6 cm/s]. A similar exercise was performed for Papp where measurements were transformed using log10(Papp) to scatter the response variable. The Papp dataset was also randomly split into a training set (70%, 3087 compounds) and a test set (30%, 1324 compounds). There was no bias in the separation of the dataset (p > 0.05; Wilcoxon signed-rank test), and the distributions are shown in Figure 2.

Figure 3. Distribution of Papp dataset

The x-axes in the figures above (left to right) show Papp values, converted Papp, molecular weight, and SlogP. The y-axis indicates the number of samples. The blue and orange bars represent the training and test sets distribution, respectively.

2.5. Descriptor preparation

We performed data preparation steps for descriptors of compounds in the training set. First, descriptors with nan were removed. Next, descriptors with small variance were removed using the VarianceThreshold function (threshold=1.0), and the retained descriptors were normalized using the StandardScaler function. Finally, effective descriptors for prediction were selected using ElasticNet. ElasticNetCV (cv=10, max_iter=1000000) selected proper alpha (penalty) and l1_ratio (norm ratio), using which we performed ElasticNet to retain proper descriptors. The number of retained descriptors and decided parameters for each dataset are shown in Tables 1 and 2. These descriptors were also used for test compounds.

Table 1. Summary of descriptor selection for Fa prediction

Dataset 1D and 2D 1D, 2D and 3D
Fa Number of descriptor 199 197
Parameter: alpha 0.10207 0.12048
Parameter: l1_ratio 0.1 0.1

Table 2. Summary of descriptor selection for Papp prediction

Dataset 1D and 2D 1D, 2D and 3D
Papp Number of descriptor 243 277
Parameter: alpha 0.0006 0.0006
Parameter: l1_ratio 1.0 0.99

2.6. Model construction

To determine the suitable algorithm for two parameter predictions, we compared Multi Linear Regression (MLR), Support Vector Regression (SVR), and Random Forest Regression (RFR). MLR is a simple regression model. SVR is a kernel function model where we employed the frequently used radial basis function kernel. Here, three parameters need to be optimized, which are C: regularization parameter, epsilon: insensitive zone, and gamma: kernel coefficient. RFR is an ensemble model based on the decision tree method that requires us to optimize the following parameters: n_estimators: the number of trees, max_depth: the maximum depth, min_samples_split: minimum number of samples required to split an internal node. The suitable parameters were selected using GridsearchCV (cv=10, score=’neg-mean-squared-error') and are shown in Tables 3 and 4.

Table 3. Summary of selected parameters for Fa prediction

Dataset 1D and 2D 1D, 2D and 3D
Fa SVR

C: 32

epsilon: 1

gamma: 0.00195

C: 32

epsilon: 2

gamma: 0.00098

RFR

n_estimators: 1000

max_depth: 25

min_samples_split: 2

n_estimators: 1200

max_depth: 25

min_samples_split: 2

Table 4. Summary of selected parameters for Papp prediction

Dataset 1D and 2D 1D, 2D and 3D
Papp SVR

C: 8

epsilon: 0.25

gamma: 0.00781

C: 8

epsilon: 0.125

gamma: 0.00781

RFR

n_estimators: 1500

max_depth: 65

min_samples_split: 2

n_estimators: 1750

max_depth: 45

min_samples_split: 2

2.7. Performance evaluation

The R2 is an inadequate score for nonlinear models [14]. Thus, we used the root mean squared error to compare the predictive performance of the three algorithms.

3. Results

3.1. Predictive accuracy

The constructed models were used to predict Fa and Papp in the test sets without 3D and with 3D descriptors (Tables 5 and 6). The RFR and SVR exhibited excellent performance for Fa and Papp prediction model, respectively. Furthermore, the 3D descriptors were adequate to predict Fa (Figure 3).

Table 5. Summary of Fa prediction results of root mean squared error for the test set

Dataset 1D and 2D 1D, 2D and 3D
Fa MLR 5.2702 4.0395
SVR 2.0109 2.0214
RFR 1.9103 1.9258

Table 6. Summary of Papp prediction results of root mean squared error for the test set

Dataset 1D and 2D 1D, 2D and 3D
Papp MLR 1.6241 1.6150
SVR 1.2978 1.2968
RFR 1.3541 1.3548

Figure 4. Left scatter plot shows the result of the observed and predicted Fa

The predicted Fa was calculated with RFR model constructed using 1D and 2D descriptor. The right scatter plot shows the result of the observed and predicted Papp. The predicted Papp was calculated with SVR model constructed using 1D, 2D, and 3D descriptors. The x-axis and y-axis show the values of converted Fa and Papp, respectively. The correlation coefficient scores between the observed and predicted values were 0.6861 (Fa) and 0.7019 (Papp).

3.2. Contribution of descriptors for Fa prediction

Although the generated relationship between the descriptors and predicted results were effective for compound optimization, the constructed models were complex, and it was difficult to elucidate their relationships. Shapley additive explanations (SHAP) is a useful tool to overcome this hurdle [15], where the method calculates an important value of each feature for a prediction based on game theory. In Fa prediction by RFR, low values of TopoPSA(NO) and VSA_EState are likely to increase output value (Figure 4). In contrast, the low values of SLogP and ATSC3p tend to decrease results. The calculation performed by the method was time- and memory-intensive due to which we did not compute the SHAP in Papp dataset. However, the script to calculate SHAP value is shown in supplemental information for Papp prediction.

Figure 5. Contribution of descriptors for Fa predicted by RFR

The plot shows important SHAP values for Fa prediction. The x-axis and y-axis show the values of SHAP value and descriptors, respectively. The red points represent high feature value and the blue points show low value.

4. Conclusion

Although numerous ML packages have facilitated the construction of in silico models, this ease of model construction is beyond the reach of researchers who lack expertise in programming. These wet-lab researchers are forced to rely on package documentation and websites for model construction, appropriate data preparation, and evaluation. We prepared and shared the model construction scripts to predict Fa and Papp of chemical compounds. These scripts require Jupyter environment and SDF with a 3D structure for Papp prediction. We are confident that the application of this model will accelerate the effective search for potential drug candidates and sustain drug discovery. These models, along with the protocol, have been shared in the supporting information for the reader.

Supplemental information

The lists of descriptors, scripts for Fa and Papp prediction, and other information are available in as follows:

  •    CBIJ_Esaki_et_al_Descriptor_Mordred_1D2D.csv: List of 1D and 2D descriptor
  •    CBIJ_Esaki_et_al_Descriptor_Mordred_1D2D3D.csv: List of 1D, 2D and 3D descriptor
  •    CBIJ_Esaki_et_al_ModRep_DescriptorCalculation.ipynb: Script for descriptor calculation
  •    CBIJ_Esaki_et_al_ModRep_FaModel.ipynb: Script for Fa model construction
  •    CBIJ_Esaki_et_al_ModRep_PappModel.ipynb: Script for Papp model construction
  •    CBIJ_Esaki_et_al_ModRep_SI.pdf: Protocol for model construction, summary of tuned parameters, and scatter plots of the observed and the predicted values.

Acknowledgments

This work was supported by JSPS KAKENHI (Grant Number 20K16075). We thank Editage (www.editage.jp) for English language editing.

Conflict of Interest

The authors declare no conflict of interest.

References
 
International (CC BY 4.0) : The images, videos or other third party material in this article are also included in the article’s Creative Commons license.To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/
feedback
Top