Cheaper faster drug development validated by the repositioning of drugs against neglected tropical diseases

There is an urgent need to make drug discovery cheaper and faster. This will enable the development of treatments for diseases currently neglected for economic reasons, such as tropical and orphan diseases, and generally increase the supply of new drugs. Here, we report the Robot Scientist ‘Eve’ designed to make drug discovery more economical. A Robot Scientist is a laboratory automation system that uses artificial intelligence (AI) techniques to discover scientific knowledge through cycles of experimentation. Eve integrates and automates library-screening, hit-confirmation, and lead generation through cycles of quantitative structure activity relationship learning and testing. Using econometric modelling we demonstrate that the use of AI to select compounds economically outperforms standard drug screening. For further efficiency Eve uses a standardized form of assay to compute Boolean functions of compound properties. These assays can be quickly and cheaply engineered using synthetic biology, enabling more targets to be assayed for a given budget. Eve has repositioned several drugs against specific targets in parasites that cause tropical diseases. One validated discovery is that the anti-cancer compound TNP-470 is a potent inhibitor of dihydrofolate reductase from the malaria-causing parasite Plasmodium vivax.


Yeast strains and plasmid constructs
The list of strains constructed is given in Table S1. These are described in Bilsland et al. (1). Briefly, fluorescent plasmids were constructed by replacing the coding region of yEmRFP from yEpGAP-Cherry (2) with Venus or Sapphire (3) and replacing the URA3 marker with LEU2. The strain expressing the drug-resistant P.vivax DHFR ( PvR dhfr) was constructed by mutating the following sites of the target enzyme: S58R, S117N and I173L. The plasmid was transformed into a yeast strain with a dfr1∆/DFR1 pdr5∆/PDR5 BY4743 background. The strain was sporulated and MAT haploids were selected for drug screens. The drug resistant PfR dhfr strain is a triple mutant for residues N51I, C59R & S108N, the required changes to the wildtype sequence having been made by site-directed mutagenesis.
These three fluorescent proteins have very distinct excitation and emission patterns: mCherry (ex 580nm; em 612nm), Sapphire (ex 405nm; em 510 nm) and Venus (ex 500nm; em 540 nm).  Table S1. Yeast strains and plasmids The Yeast strains were cultured in YNB-glucose (0.68% yeast nitrogen base without amino acids, 2% ammonium sulphate, 2% glucose) with the relevant supplements for all assays. The initial plan of Eve was given in (4) (Fig. S1), and hardware and lov level software are described in the main text.

Establishment of the mixed cultures
Pre-cultures were grown to stationary phase and 1mL of each culture was inoculated into 100mL of fresh medium. Pools were incubated at 30°C, with shaking, for 4h to ensure exponential growth. Doxycycline (5µg/mL) was then added to the culture to reduce expression of the target enzyme (1). The culture was attached to a Thermo Combi multidrop within the Eve work cell. The culture was stirred continuously and maintained at 23°C during assay plate set up.

Screening and hit confirmation, using the Robot Scientist Eve
Screens were performed by the Robot Scientist Eve using mixed cultures described above and either the ~1,600 FDA-and foreign-approved drugs from the Johns Hopkins University Clinical Compound Library or the 14,400 compound Maybridge Hitfinder library. Strains were grown in competition in the presence of a library compound, as discussed above, and the relative growth rates used to estimate the activity of the drug against the parasite target.
A typical Maybridge library screen consists of a set of 45 384-well plates, each well containing a pool of three yeast strains harboring either a parasite drug target or its human counterpart. Different Maybridge compounds were added to each of 320 wells on each plate. 20 of the remaining wells were used for positive controls (pyrimethamine, trimethoprim, methotrexate, raltitrexed monohydrate, pemetrexed disodium), and 44 for negative controls, with strains growing in an equivalent concentration of DMSO. Each well was inoculated with 50 L of the pooled yeast culture (final compound concentration of 10 M). Plates were incubated for 40 hours, and fluorescence measurements taken every 90 minutes.

Growth curve parameter derivation
Growth curves were fit to the time course, and growth parameters A-P (Fig. S2) derived (5,6). Figure S2. Typical growth curve, and growth parameters extracted from the curve.

Reproducibility of growth curve parameters
The variation of doubling time (DT) amongst the negative controls within a plate was used to estimate in-plate repeatability, and between-plate reproducibility was estimated using DT of both the negative and positive controls. Example data from screen TS3 (Hs, Pv, PfR) is shown in Table S2  DT was then used to identify compounds causing significantly reduced growth (> two standard deviations) relative to the negative controls. For screen TS3, approximately 600 compounds were identified as potential hits in this way. These were then classified according to the remaining growth parameters (Fig. S2), as (a) autofluorescent; (b) cross-inhibited and (c) hit compounds (Fig. S3) using the following criteria: (a) Autofluorescent compounds -fluorescence of the drug itself interfering with the growth assay 1. High fluorescence at one or more wavelengths throughout the run 2. High initial fluorescence (>8% higher than the negative control value) at one or more wavelengths (b) Compounds with evidence of cross-inhibition 1. Low intensity fluorescence for all strains 2. Low growth for all strains 3. Long lag time for all strains Note: Compounds with these properties can either be toxic to yeast (in that case, they would appear as toxic in all of the screens, irrespective of the target), problem wells (technical/experimental error), or compounds that inhibit yeast growth when either the parasite enzyme or its human ortholog are expressed.
(c) Hit compounds -candidate anti-parasitic drugs 1. Reduced yield for a parasite strain relative to the human strain 2. Reduced growth rate for a parasite strain relative to the human strain  Figure S3. Typical forms of growth curves observed in Eve drug screens, and the classification of compounds based on the curve properties.
Derivation of Rules for deciding library-screening hits 325 of these initial candidate hit compounds, based on abnormal DT, were categorised manually, resulting in:  83 possibly toxic compounds  16 autofluorescent compounds  57 strong and 64 weak hits against PvDHFR  9 weak hits against PfRDHFR  67 inactive compounds Based on this primary data, Weka 3.6.2, C4,5 (J48) decision trees (7,8) were used to determine an optimal set of rules for automatically classifying growth curves.
The rules defined by the decision trees were based on a subset of the growth parameters (Fig. S2). These were: the strain 'yield ratio' -the change in fluorescence when grown in the presence of a drug, relative to the negative control; the initial strain fluorescence (relative to the change in fluorescence of the negative control); doubling time relative to the negative control; and lagtime2 relative to the negative control. The combined growth of all strains within a well was also employed in the decision trees to distinguish cross-inhibited compounds.
The decision tree process produced very similar rules when by considering each of the Hs, Pv and PfR datasets individually. These were combined into a generalised set of rules which were used in all subsequent screens (Table S3). Table S3. Machine learning rules derived to identify potential hit compounds with the Eve mass screen.

Analysis and confirmation of intelligent screening data
After performing a mass screen and classifying compounds according to the rules defined above, confirmation screens were performed on the potential hit compounds. The same protocols were employed during the intelligent screening process. Machine learning rules for the cherry-picked screens were developed using Weka (J48 decision trees) as above. Rules were validated by comparing machine-learning classifications with inspection of the growth curves, which resulted in near-complete agreement for all screens. If more than four of the 40 individual cherry curves are hits, then the compound is active against the cherry-labelled (Hs) target, and is classed as possibly cross-inhibited.

Machine Learning
To form QSAR hypotheses Eve uses a Gaussian process model (9), with the molecules represented using binary fingerprints representing all linear paths of up to seven atoms, computed using Open Babel (10). The inputs for learning were OpenBabel FP2 fingerprints (using the 0/7 configuration) for compound SMILES codes (training set and unknowns), and yield ratio differences between target and human strains, e.g (strain_yield_ratio HsDHFR -strain_yield_ratio PvDHFR ).

Implementing Active Learning loops
The data analysis and subsequent rules were explained above and were combined with the following to enable an autonomous process:  a method to select 'hit' compounds (either by active learning or by selection for a confirmation screen based on the decision tree rules)  a method to combine the confirmation/intelligent screen assays complex cherrypick results with the simpler mass screen data, to build the information source for the active learning feedback loop (Fig. S4). Figure S4. Eve's strategy for improving the efficiency of library screening by the application of machine learning.
The first full experimental tests of the active learning loop were conducted by splitting the screen data set for TS6 (comprising the heterologous DHFR yeast strains for Plasmodium falciparum and P.vivax, and that of humans), and using 4800 compounds as a training set. The yield ratios of the HsDHFR and PvDHFR strains were passed to the selection algorithm, together with fingerprints of the remaining 9600 compounds. The results from the first cherry-picking round (n=96; 12 plates of 8 compounds per plate; 8 replicates of 6 concentrations) were then added to the original data set, and a second cherry-picking round conducted. To examine different approaches to the problem of combining cherry-picking and mass screening data, seven versions of weighted cherry-picking data (strain_yield_ratio HsDHFR -strain_yield_ratio PvDHFR ) were tested in the second round: i. All replicates for each compound ii. All replicates multiplied by 10/conc iii. All replicates multiplied by log(10/conc) iv.
Mean of replicates at each concentration for each compound v. Mean of replicates multiplied by 10/conc vi.
No additional data (i.e. the next best 96 compounds under loop1 conditions) When the ML rules were applied to TS6, 282/14099 compounds (2.0%) were identified as active against the Pv strain. The curves for the cherry-picking data were empirically categorised, and evaluated against the weighting options, with the log-weighted option (vi) being deemed optimal. Table S5. Methods for integrating cherry-picking screen data into Eve's active learning algorithm.

Active Learning simulations
The active learning loop was run through three iterations; an initial set of 4800 compounds was screened (single iteration, 10 µM), and three loops of 96 cherrypicked compounds (8 replicates, at a range of concentrations) were selected. The mean log.-weighted cherry-picking data is cycled back into the training set. In addition to physical active-learning screens, simulations were performed, using the empirical data for each compound in the whole-library screen dataset to simulate logweighted cherry-picking scores for the compound. Experiment and simulation are compared in Table S6.

Econometric modelling
The model To determine the utility of Eve for drug discoveryi.e. the range of conditions for which using a Robot Scientist to guide candidate compound selection is economically advantageous compared with performing a standard whole-library screenwe developed the econometric model presented in Figure 2 in the main text (and below).  The net utility is made up of three components: the cost saving due to not screening Nm compounds which, based on the QSAR learning, are unlikely to be hits, minus the opportunity cost due to not finding any hits (Nx) that might be present in this unscreened set (Uh >> Tc + Cc ), minus the cumulative cost of the number of active learning cycles performed [the cost of cherry-picking Ne compounds, where (Tc + Cc ) > (Tm + Cm )].

Active k-optimisation strategy
One of the major goals for Eve's data analysis was to build algorithms to predict active compounds. The array of information contained in assay data includes raw data for yeast target growth profiles, labelled classifications for activity, toxicity etc., and structural representations of the compounds.
The prototype method for this work is based on the active k-optimisation strategy (11). This strategy is introduced for machine learning as a way of finding and ranking the k best alternatives for evaluation, using Gaussian process to provide a mechanism for developing this model. The process finds more than one target (other than the optimal solution) to take into the next step of HTE. The work is based on having a finite library of examples, of which the results from the known ones can be used to pick the best unknowns for evaluation.
In this particular approach, the goal is to pick targets that have the best chance of success (the maximum predicted strategy) for comparison to several other existing approaches (6,12,13). The lower confidence bound criterion (optimistic) (14), selecting the sample with the highest probability of improving the current solution (most probable improvement, MPI) and efficient global optimisation (EGO) to give maximum expected improvement.
Specific application of the active k-optimisation strategy to the drug screening process is provided in (15); this describes the analysis of the NCI60 dataset (US National Cancer Institute, 60 anticancer drug screen). The full techniques upon which this strategy is based are given in (9).

Cherry-picking simulations using active k-optimisation
The inputs for the prototype AL method are OpenBabel FP2 fingerprints (using the 0/7 configuration) (16,17) for compound SMILES codes (18,19) (training set and unknowns), and simple growth differences between target and human strains, i.e.
The method is designed to identify the next best compound to test from a given library, whilst avoiding selection of a compound very similar to any in the training set; it is expected to balance exploration of the full chemical space of the library, with exploitation of areas most likely to contain active examples.

Econometric modelling using active k-optimisation simulation data
The active k-optimisation AL algorithm was applied to the seed input data and the unknown compound SMILES codes; simulated learning curves were produced for each parasite strain using the proxy confirmation data (see section 4.2). The progression of these learning curves was then compared to the base case of a linear progression throughout the screen in accordance to the utility equation (fig. S5). For each 96-compound loop, the number of proxy confirmed hits and compounds screened to date (Ne) were applied to the utility equation, together with the fixed utility and cost terms. An example of the resultant 2D plot for the PvDHFR strain is shown in fig. S6.  0  n+1  226  244  260  279  291  n+2  228  247  263  280  n+3  229  250  264  281  n+4  231  251  268  283  n+5  233  252  271  284  n+6  235  253  273  285  n+7  238  254  274  286  n+8  240  255  275  288  n+9  242  256  277  289  n+10 243 259 278 290 Figure S6: Hits found in TS3 PvDHFR simulation (red), base case (black), resultant econometric utility model (blue).
The cost of a cycle of ML includes both the time and cost of the computing power, and the cost of testing a 96-compound batch in the cherry-picking assay. Based on the active k-optimisation simulated cherry-picking cycles, these terms were also calculated for examples of each parasite target (figure S7).
A utility landscape was also constructed for the TS3 PvDHFR parasite data across a range of ML efficiencies and drug economic values (figure 2c); the value of a hit compound for this study ranged from $3K to $25K, based on a broad estimate of the number of hits required to give sufficient drug-like lead compounds to commence lead optimisation studies. Variation in the time-cost ratio comparing mass to intelligent screening (Tc/Tm) was studied, and utility versus cost of compound loss during cherry-picking (Uh/Cc) was also evaluated.

Ontological Description of the Data
Tables S8 and S9 contain the list of the relations, classes, and their URIs used for the semantic model shown at Fig. S8 has-input RO An entity X participates in a process Y and X is present at the beginning of the process Y.
hasoutput RO An entity X participates in a process Y and X is present at the end of the process Y.

fluorophore PPI
A fluorophore is a component of a molecule which causes a molecule to be fluorescent. It is a functional group in a molecule which will absorb energy of a specific wavelength and re-emit energy at a different (but equally specific) wavelength. The amount and wavelength of the emitted energy depend on both the fluorophore and the chemical environment of the fluorophore.

In vitro Enzyme Assays
Materials and Methods TNP-470 was purchased from Sigma-Aldrich.

Purification of DHFRs
Fresh overnight culture from a single colony of E. coli BL21(DE3)pLysS cells harboring wild-type PvDHFR (20) and human DHFR were used to inoculate 1-4 liters of LB medium supplemented with 100 g/ml ampicillin. The culture was grown at 37˚C until OD 600 reached 0.5-0.6, at which time the expressions of PvDHFR and human DHFR were initiated by addition of IPTG at a final concentration of 0.4 mM. The culture was allowed to grow with vigorous shaking at 20˚C for additional 18-20 h prior to harvesting by centrifugation at 10,000 g for 10 min at 4˚C. The cell pellet was washed once with ~200-250 mL of cold phosphate-buffered saline, followed by resuspending in ~20-50 mL of buffer A (20 mM potassium phosphate buffer, pH 7.0, 0.1 mM EDTA, 10 mM DTT). The cells were then disrupted by two cycles of French press at 18,000 psi. The clear extract obtained after centrifugation at 30,000 g for 1 h at 4˚C was applied onto an MTX-Sepharose CL-6B affinity column (1.5 x 5.0 cm), and the enzymes were affinity purified according to the procedure described previously (21).

Enzyme assay and Inhibition Studies
The activity of DHFR was determined spectrophotometrically by monitoring the rate of decrease in absorbance at 340 nm (22,23). The standard DHFR assay (1 ml) performed in 1-cm path-length cuvette contained 1x DHFR buffer (50 mM TES, pH 7.0, 75 mM β-mercaptoethanol, 1 mg/mL bovine serum albumin), 100 M H 2 folate, 100 M NADPH, and 0.01 units of affinity-purified enzyme. The reaction was initiated with H 2 folate. One unit of DHFR activity is defined as the amount of enzyme that produces 1 mole of product per minute at 25 C. Inhibition studies were carried out by including increasing amounts of the tested compounds in the assay reaction. The reaction was initiated with H 2 folate, and the decrease in absorbance at 340 nm was monitored as mentioned above.