Data
Data for the SARS-CoV-2 cytopathic effect (CPE), the Vero E6 (host cell) cytotoxicity counterscreen, the Spike-ACE2 protein-protein interaction (AlphaLISA), TruHit (AlphaLISA counterscreen), ACE2 enzymatic inhibition conterscreen, and 3CL enzymatic inhibition were obtained from the NCATS COVID19 portal.4,5 Assay data were mapped to DrugCentral2020 to retrieve DrugCentral IDs, SMILES strings, and drug names. The experimental results from NCATS were mined to define activity and significance classes. The procedure for data extraction, wrangling, and post-processing are detailed as follows.
Mining NCATS COVID19 Data
All operations described below were performed using the Knime35 platform. NCATS data associated with the aforementioned assays was downloaded from the COVID19 portal.4,5 The files contained over 23,000 compounds generated by high-throughput screening (HTS) experiments. When possible, each compound was cross-linked to drugs annotated in DrugCentral, to separate approved drugs from the rest of the tested compounds. Matching NCATS compounds to DrugCentral was conducted in three sequential steps: by InChI (International Chemical Identifier),36 by synonym (name), and by matching CAS (Chemical Abstracts Service) Registry Numbers. First, NCATS molecules were matched by InChI. Molecules that did not match were then queried by drug name and associated synonyms, as annotated in DrugCentral. Finally, if not matched by either InChI or name, molecules were matched by CAS number. If none of the above steps resulted in a match, then the molecule in question was not classified as an approved drug. At the end of this process, 4,954 unique molecules (2,273 approved drugs and 2,681 chemicals) were stored. Whenever possible, for approved drugs and other chemicals, SMILES were retrieved from DrugCentral and from PubChem, respectively. Otherwise, the original SMILES strings were retained. Bioactivity data was mined according to the “curve class” and “maximum response” parameters.37 The “curve class” evaluates the quality of the dose-response curve and, consequently, the quality of the measurement. The “maximum response” is the maximum response value detected during the experiment. The “ACTIVITY CLASS” and a “SIGNIFICANCE CLASS” were defined using criteria described in Table 4 and Table 5, respectively. As a final data wrangling step, all compounds were categorized, and assay data grouped to have a unique record per molecule for each assay. When more than one assay was measured for the same molecule, only the datapoint with the best curve class was retained.
Table 4: Criteria applied to each experimental dataset that define the various activity categories. Depending on the assay type, MAX RESPONSE could have either negative or positive values. In any case, since a greater (absolute) value always represents high activity, the absolute value of MAX RESPONSE is used in the definition of activity classes.
Cutoff
|
ACTIVITY CLASS
|
MAX RESPONSE > 66
|
HIGH
|
33 ≥ MAX RESPONSE ≤ 66
|
MODERATE
|
MAX RESPONSE < 33
|
LOW
|
Table 5: The significance classes represent the quality of the measurement and are defined in accordance with the different shape of the dose-response curve (CURVE CLASS2).37 (a) A type 4 curve is a flat curve, which indicates that the compound does not show any activity.
CURVE CLASS2
|
SIGNIFICANCE CLASS
|
1.1 1.2 2.1 2.2
|
HIGH
|
1.3 1.4 2.3 3.4
|
MODERATE
|
3
|
LOW
|
4
|
INACTIVEa
|
5
|
INCONCLUSIVE
|
Data Filtering
For each assay, the data was labeled as positive and negative. The compounds with “LOW” activity class were treated as negative, whereas “HIGH” and “MODERATE” were treated as positive compounds. Molecular structures were standardized using five different filters using RDKit38 features: (i) SMILES strings were converted into canonical SMILES strings. Some of the SMILES were not converted SMILES, were discarded. (ii) RDKit Salt Stripper was implemented to obtain the salt stripped molecules. The “donotRemoveEverything” feature helps by leaving the last salt when the entire canonical SMILES string contains only the salts, (iii) The RDKit “Uncharger” feature neutralizes the molecules by adding or removing hydrogen atoms, (iv) The canonical tautomer was used obtained from RDKit, and (v) logP and logS filters were defined as follows: (logP < 1 and logP > 9) and (logS > -3 and logS < -7.5). Such thresholds were arbitrarily defined to maximize the number of negative data excluded while minimizing the positive data excluded. However, for the 3CL dataset, logP and logS filters were not implemented. The final task was the removal of duplicate compounds resulting from desalting, neutralizing, and standardizing tautomers.
Molecular Descriptors
A total of 22 features of three distinct types (fingerprints-based, pharmacophore-based, and physicochemical descriptors based) were implemented. Fingerprints were converted into a bit vector of either 1,024 or 16,384 lengths. Pharmacophore type was also a bit vector of size 2,692, whereas RDKit and VolSurf+ descriptors were of length 200 and 128, respectively.
Fingerprints-based: This includes the circular, path-based, and substructure keys.39,40 Circular fingerprints include the extended-connectivity fingerprints (ECFPx) and feature-connectivity fingerprints (FCFPx), where x is 0, 2, 4, and 6 are the bond length or diameter for each circular atom environment. ECFP consists of the element, number of heavy atoms, isotope, number of hydrogen atoms, and ring information, whereas FCFP consists of pharmacophore features.
Substructure keys Avalon and the public Molecular ACCess System (MACCS) are two distinct types of fingerprints that are substructure keys. The Avalon fingerprint, used here, is a bit vector of size 1,024. It includes feature classes such as atom count, atom symbol path, augmented atom, and augmented symbol path, etc. MACCS structural keys are 166-bit structural key descriptors. Each bit here is associated with a SMARTS pattern and belongs to the dictionary-based fingerprint class. Path-based fingerprints include RDKx (where x is 5, 6, 7), topological torsion (TT), HashTT, atom pair (AP), and HashAP. The size of each fingerprint is 1024 for all of them.
A longer version of the fingerprint, of 16384 bits, was also used for comparison. This longer version is represented by the prefix” L”: LAvalon, LECFP6, LECFP4, LFCFP6 and LFCFP4.
Pharmacophore-based: Topological pharmacophore atomic triplets fingerprints (TPATF) were obtained using maychemtools.41 TPATF describes the ligand sites that are necessary for molecular recognition of a macromolecule or a ligand, and passes that information to the ML model to be trained. Ligand SMILES strings were passed through a Perl script to generate TPATF. The basis sets of atomic triplets were generated using two different constraints (i) triangle rule, i.e., the length of each side of a triangle cannot exceed the sum of the lengths of the other two sides; and (ii) elimination of redundant pharmacophores related by symmetry. The default pharmacophore atomic types Hydrogen Bond Donor (HBD), hydrogen bond acceptor (HBA), positively ionizable (PI), negatively ionizable (NI), H (hydrophobic), and Ar (aromatic) were used during generation of TPATF.42
Physicochemical descriptors: This includes the RDKit molecular descriptors and VolSurf+ descriptors. For RDKit descriptors, a set of 200 descriptors were used, which were obtained from RDKit.38 They are either experimental properties or theoretical descriptors, which are e.g. molar refractivity, logP, heavy atom counts, bond counts, molecular weight, topological polar surface area.
A total of 128 descriptors were obtained using VolSurf+ software. VolSurf+ is a computational approach aimed at describing the structural, physicochemical and pharmacokinetic features of a molecule starting from a 3D map of the interaction energies between the molecule and chemical probes (GRID-based molecular interaction fields, or MIFs).43 VolSurf+ compresses the information present in MIFs into numerical descriptors, which are simple to use and interpret.16
Model Development
For each classifier, several machine learning models were developed, employing three categories of features and 22 distinct machine learning algorithms from the scikit-learn package.18 Figure 8 briefly shows the workflow of the model generation. The three different categories of features employed were 1) chemical fingerprints, 2) physicochemical descriptors, and 3) topological pharmacophore descriptors. For fingerprint-based descriptors, 19 different RDKit fingerprints were tested. For physicochemical descriptors, Volsurf+ and RDKit descriptors were employed. For Topological pharmacophore descriptors, TPATF fingerprints from Mayachemtools were employed. For each model, input data was split into a 70% training set, 15% validation set, and 15% test set using a stratified sampling. Table 6 depicts the number of compounds used in training, validation, and test sets for each model. Compounds discussed in the External Predictivity section (Figure 7) were not part of this workflow and are not included in Table 6.
All of the models were trained leveraging 22 different classifiers available in scikit-learn. Initially, 968 different models were trained using 22 classifiers (Complement Naive Bayes, Extreme Gradient Boosting, KNeighbors, Gradient Boosting, Perceptron, One Vs Rest , Extra-Tree, Ridge, One Vs One, Bagging, Random Forest, Output Code, Passive Aggressive, Linear SVC, Stochastic Gradient Descent, Logistic Regression, Extra Trees, Multinomial Naive Bayes, Ada Boost, Decision Tree, Nearest Centroid, Multi-layer perceptron) and 22 distinct features (ECFP0, ECFP2, ECFP4, LECFP4, ECFP6, LECFP6, FCFP2, FCFP4, LFCFP4, FCFP6, LFCFP6, RDK5, RDK6, RDK7, Avalon, LAvalon, MACCS, HashTT, HashAP, VolSurf+, TPATF, and RDKit descriptors) based on the default configurations of the classifiers. Finally, the best suited combination of classifiers and features were selected for hyperparameter tuning.
Table 6: Summary of datasets used for each model.
Type
|
Count
|
Total
Actives
|
Training
count
|
Validation
count
|
Test
count
|
Training
Actives
|
Validation
Actives
|
Test
Actives
|
CPE
|
736
|
368
|
532
|
102
|
102
|
266
|
51
|
51
|
cytotox
|
1662
|
877
|
1163
|
249
|
250
|
581
|
125
|
125
|
AlphaLisa
|
1260
|
630
|
882
|
189
|
189
|
441
|
94
|
95
|
TruHit
|
1658
|
829
|
1161
|
248
|
249
|
580
|
124
|
125
|
ACE2
|
206
|
103
|
144
|
31
|
31
|
72
|
15
|
16
|
3CL
|
442
|
221
|
309
|
66
|
67
|
154
|
33
|
34
|
Similarity Search
We used an ECFP4 bit vector fingerprint with 1024 bits, and TC calculations, for the fingerprints present in the database along with that of a query molecule, are computed on the fly. The Tanimoto coefficient represents the overlap of features between molecules as the ratio of the number of common features to the total number of features in each fingerprint. The coefficient ranges from 0 to 1, with 1 corresponding to identical fingerprints.
Hyperparameter Optimization
The best performing default models based on validation sets were selected for hyperparameter optimization. These models were optimized using a grid search method. Scikit-learn provides the package for grid search hyperparameter optimization using cross validation. However, it is slow and does not offer a grid search with a separate validation set. Thus, we used a freely available software “hypopt” for the hyperparameter tuning of each model.44 The models were optimized and selected based on the validation F1 score. The outperforming models were saved and used for the evaluation of external datasets.
Evaluation Metrics
The model performances were evaluated by five distinct evaluation metrics available in scikit-learn: accuracy, recall, precision, F1 score, and area under the receiver operating characteristic (ROC) curve (AUC). In general, the accuracy measures the total number of correct predictions among the total numbers of instances evaluated. The recall was used to measure the fraction of total positives that are correctly classified, whereas precision estimates the fraction of total positives that are correctly predicted from the total predicted positives. F1 score is the harmonic mean of recall and precision. The score ranges from 0 to 1 for each of the metrics, where 1 is a perfect score. Moreover, AUC was computed for each of the classifiers. A perfect classifier yields an AUC of 1, whereas a random model has an expected AUC of 0.5.
Implementation and Accessibility
Web Portal
REDIAL-2020 is available online at http://drugcentral.org/Redial.
Code Availability
All of the codes and the trained models are available at https://github.com/sirimullalab/redial-2020.