Quantitative Systems Pharmacology and Machine Learning: A Match Made in Heaven or Hell?

As pharmaceutical development moves from early-stage in vitro experimentation to later in vivo and subsequent clinical trials, data and knowledge are acquired across multiple time and length scales, from the subcellular to whole patient cohort scale. Realizing the potential of this data for informing decision making in pharmaceutical development requires the individual and combined application of machine learning (ML) and mechanistic multiscale mathematical modeling approaches. Here we outline how these two approaches, both individually and in tandem, can be applied at different stages of the drug discovery and development pipeline to inform decision making compound development. The importance of discerning between knowledge and data are highlighted in informing the initial use of ML or mechanistic quantitative systems pharmacology (QSP) models. We discuss the application of sensitivity and structural identifiability analyses of QSP models in informing future experimental studies to which ML may be applied, as well as how ML approaches can be used to inform mechanistic model development. Relevant literature studies are highlighted and we close by discussing caveats regarding the application of each approach in an age of constant data acquisition. SIGNIFICANCE STATEMENT: We consider when best to apply machine learning (ML) and mechanistic quantitative systems pharmacology (QSP) approaches in the context of the drug discovery and development pipeline. We discuss the importance of prior knowledge and data available for the system of interest and how this informs the individual and combined application of ML and QSP approaches at each stage of the pipeline.

Free Drug Theory - No Longer Just a Hypothesis?

The Free Drug Hypothesis is a well-established concept within the scientific lexicon pervading many areas of Drug Discovery and Development, and yet it is poorly defined by virtue of many variations appearing in the literature. Clearly, unbound drug is in dynamic equilibrium with respect to absorption, distribution, metabolism, elimination, and indeed, interaction with the desired pharmacological target. Binding interactions be they specific (e.g. high affinity) or nonspecific (e.g. lower affinity/higher capacity) are governed by the same fundamental physicochemical tenets including Hill-Langmuir Isotherms, the Law of Mass Action and Drug Receptor Theory. With this in mind, it is time to recognise a more coherent version and consider it the Free Drug Theory and a hypothesis no longer. Today, we have the experimental and modelling capabilities, pharmacological knowledge, and an improved understanding of unbound drug distribution (e.g. Kpuu) to raise the bar on our understanding and analysis of experimental data. The burden of proof should be to rule out mechanistic possibilities and/or experimental error before jumping to the conclusion that any observations contradict these fundamentals.

Identifying and characterising the impact of excitability in a mathematical model of tumour-immune interactions.

We study a five-compartment mathematical model originally proposed by Kuznetsov et al. (1994) to investigate the effect of nonlinear interactions between tumour and immune cells in the tumour microenvironment, whereby immune cells may induce tumour cell death, and tumour cells may inactivate immune cells. Exploiting a separation of timescales in the model, we use the method of matched asymptotics to derive a new two-dimensional, long-timescale, approximation of the full model, which differs from the quasi-steady-state approximation introduced by Kuznetsov et al. (1994), but is validated against numerical solutions of the full model. Through a phase-plane analysis, we show that our reduced model is excitable, a feature not traditionally associated with tumour-immune dynamics. Through a systematic parameter sensitivity analysis, we demonstrate that excitability generates complex bifurcating dynamics in the model. These are consistent with a variety of clinically observed phenomena, and suggest that excitability may underpin tumour-immune interactions. The model exhibits the three stages of immunoediting - elimination, equilibrium, and escape, via stable steady states with different tumour cell concentrations. Such heterogeneity in tumour cell numbers can stem from variability in initial conditions and/or model parameters that control the properties of the immune system and its response to the tumour. We identify different biophysical parameter targets that could be manipulated with immunotherapy in order to control tumour size, and we find that preferred strategies may differ between patients depending on the strength of their immune systems, as determined by patient-specific values of associated model parameters.

Complete Substrate Inhibition of Cytochrome P450 2C8 by AZD9496, an Oral Selective Estrogen Receptor Degrader.

AZD9496 ((E)-3-(3,5-difluoro-4-((1R,3R)-2-(2-fluoro-2-methylpropyl)-3-methyl-2,3,4,9-tetrahydro-1H-pyrido[3,4-b]indol-1-yl)phenyl)acrylic acid) is an oral selective estrogen receptor degrader currently in clinical development for treatment of estrogen receptor-positive breast cancer. In a first-in-human phase 1 study, AZD9496 exhibited dose nonlinear pharmacokinetics, the mechanistic basis of which was investigated in this study. The metabolism kinetics of AZD9496 were studied using human liver microsomes (HLMs), recombinant cytochrome P450s (rP450s), and hepatocytes. In addition, modeling approaches were used to gain further mechanistic insights. CYP2C8 was predominantly responsible for biotransformation of AZD9496 to its two main metabolites whose rate of formation with increasing AZD9496 concentrations exhibited complete substrate inhibition in HLM, rCYP2C8, and hepatocytes. Total inhibition by AZD9496 of amodiaquine N-deethylation, a specific probe of CYP2C8 activity, confirmed the completeness of this inhibition. The commonly used substrate inhibition model analogous to uncompetitive inhibition fit poorly to the data. However, using the same model but without constraints on the number of molecules occupying the inhibitory binding site (i.e., nS1ES) provided a significantly better fit (F test, P< 0.005). With the improved model, up to three AZD9496 molecules were predicted to bind the inhibitory site of CYP2C8. In contrast to previous studies showing substrate inhibition of P450s to be partial, our results demonstrate complete substrate inhibition of CYP2C8 via binding of more than one molecule of AZD9496 to the inhibitory site. As CYP2C8 appears to be the sole isoform catalyzing formation of the main metabolites, the substrate inhibition might explain the observed dose nonlinearity in the clinic at higher doses.

Structural identifiability for mathematical pharmacology: models of myelosuppression.

Structural identifiability is an often overlooked, but essential, prerequisite to the experiment design stage. The application of structural identifiability analysis to models of myelosuppression is used to demonstrate the importance of its considerations. It is shown that, under certain assumptions, these models are structurally identifiable and so drug and system specific parameters can truly be separated. Further it is shown via a meta-analysis of the literature that because of this the reported system parameter estimates for the "Friberg" or "Uppsala" model are consistent in the literature.

Oxygen-Driven Tumour Growth Model: A Pathology-Relevant Mathematical Approach.

Xenografts--as simplified animal models of cancer-differ substantially in vasculature and stromal architecture when compared to clinical tumours. This makes mathematical model-based predictions of clinical outcome challenging. Our objective is to further understand differences in tumour progression and physiology between animal models and the clinic. To achieve that, we propose a mathematical model based upon tumour pathophysiology, where oxygen--as a surrogate for endocrine delivery--is our main focus. The Oxygen-Driven Model (ODM), using oxygen diffusion equations, describes tumour growth, hypoxia and necrosis. The ODM describes two key physiological parameters. Apparent oxygen uptake rate (k'R) represents the amount of oxygen cells seem to need to proliferate. The more oxygen they appear to need, the more the oxygen transport. k'R gathers variability from the vasculature, stroma and tumour morphology. Proliferating rate (kp) deals with cell line specific factors to promote growth. The KH,KN describe the switch of hypoxia and necrosis. Retrospectively, using archived data, we looked at longitudinal tumour volume datasets for 38 xenografted cell lines and 5 patient-derived xenograft-like models. Exploration of the parameter space allows us to distinguish 2 groups of parameters. Group 1 of cell lines shows a spread in values of k'R and lower kp, indicating that tumours are poorly perfused and slow growing. Group 2 share the value of the oxygen uptake rate (k'R) and vary greatly in kp, which we interpret as having similar oxygen transport, but more tumour intrinsic variability in growth. However, the ODM has some limitations when tested in explant-like animal models, whose complex tumour-stromal morphology may not be captured in the current version of the model. Incorporation of stroma in the ODM will help explain these discrepancies. We have provided an example. The ODM is a very simple -and versatile- model suitable for the design of preclinical experiments, which can be modified and enhanced whilst maintaining confidence in its predictions.

Tumour stromal morphology impacts nanomedicine cytotoxicity in patient-derived xenografts.

It is challenging to evaluate how tumour pathophysiology influences nanomedicine therapeutic effect; however, this is a key question in drug delivery. An advanced analytical method was developed to quantify the spatial distribution of drug-induced effect in tumours with varied stromal morphologies. The analysis utilises standard immunohistochemistry images and quantifies the frequency of positive staining as a function of distance from the stroma. Two stromal morphologies - Estuary and Tumour Island - were classified in 28 tumours from a lung cancer explant model in mice treated with liposomal doxorubicin. Analysis demonstrated that Estuary-like tumours presented a highly convoluted tumour-stroma interface, with most tumour cells in close proximity to vessels; these tumours were 8.8-fold more responsive to liposomal doxorubicin than were Tumour Island-like tumours, which were nearly unresponsive to liposomal doxorubicin. SDARS analysis allows the relative treatment effect to be assessed in tumours individually, and enables investigation of nanomedicine delivery in complex tumour pathophysiologies. FROM THE CLINICAL EDITOR: Advances in nanotechnology have brought about many novel treatment modalities for cancer. Nonetheless, there is no standard evaluation technique for tumor cells' drug response. The authors here utilized patient-derived tumour xenograft (PDTX) models to have a more translatable pre-clinical evaluation platform for nanomedicine drugs. They then used advanced imaging acquisition technique to analyze tumor stromal morphology, which they named Spatial Distribution of Apoptosis Relative to Stroma (SDARS). The findings would have significant clinical impact as it would help predict the eventual clinical drug response.

Determination of the safety and efficacy of therapeutic neutralization of tumor necrosis factor-α (TNF-α) using AZD9773, an anti-TNF-α immune Fab, in murine CLP sepsis.

OBJECTIVE AND DESIGN: TNF-α neutralization is associated with increased mortality in mouse cecal ligation puncture (CLP) models. AZD9773 is an ovine polyclonal human TNF-α immune Fab, with pharmacological properties that differ from previously studied anti-TNF-α agents. We explored the safety and efficacy of therapeutically administered AZD9773 in mouse CLP sepsis. METHODS: A moderate/severe-grade CLP model resulting in 20-30 % 5-day survival and a mild-grade CLP model resulting in ~70 % 5-day survival were established in human TNF-α transgene/murine TNF null (Tg1278/-/-) mice. TREATMENT: Mice received saline resuscitation and imipenem administration every 12 h (0-72 h post-CLP). AZD9773 (or DigiFab control) was dosed 24, 36, 48 and 60 h post-CLP. RESULTS: Therapeutic dosing of AZD9773 in moderate/severe-grade CLP resulted in significantly increased survival (>70 %) compared with DigiFab (27 %, P < 0.05). Therapeutic dosing of AZD9773 in mild-grade CLP did not significantly affect survival outcome compared with DigiFab or imipenem alone (~60-70 % survival). CONCLUSIONS: These data demonstrate that TNF-α neutralization can improve survival in moderate/severe CLP sepsis. TNF-α suppression in mild-grade models was not associated with survival benefit and did not increase 5-day mortality. These findings suggest that therapeutic benefit following TNF-α attenuation in models of sepsis may depend on model severity.

Experimental and mathematical analysis of in vitro Pitavastatin hepatic uptake across species.

1. To investigate the non-linear kinetics of in vitro hepatocyte uptake across species, the OATP substrate Pitavastatin was used as a probe. 2. Experiments were conducted at AstraZeneca (Alderley Park, Macclesfield) using freshly isolated rat, dog and human hepatocytes, utilising the "oil spin" methodology described by Hassen et al. (1996). Very few mechanistic models have previously been used to characterise the uptake process. 3. Here two candidate pharmacokinetic non-linear compartmental models are proposed. Both models have been shown to be structurally identifiable and distinghishable previously, which establishes that all unknown parameters could be identified from the experimental observations available and that input/output relationships for both the candidate models were structurally different. 4. A kinetic modelling software package, FACSIMILE (MCPA Software, Faringdon, UK), was used to obtain numerical solutions for the system equations and for parameter estimation. Model fits gave good agreement with the in vitro data and suggest the current widely accepted assumption that the rate of diffusion across the hepatocyte cell membrane is the same at both 4 °C and 37 °C is not valid, at least for Pitavastatin. Although this finding has already been proposed, this is the first time it is comprehensively debunked using statistical testing.

The design and analysis of parallel experiments to produce structurally identifiable models.

Pharmacokinetic analysis in humans using compartmental models is restricted with respect to the estimation of parameter values. This is because the experimenter usually is only able to apply inputs and observations in a very small number of compartments in the system. This has implications for the structural identifiability of such systems and consequently limits the complexity and mechanistic relevance of the models that may be applied to such experiments. A number of strategies are presented whereby models are rendered globally identifiable by considering a series of experiments in parallel. Examples are taken from the pharmacokinetic literature and analysed using this parallel experiment methodology. It is concluded that considering a series of pharmacokinetic experiments where some, but not all, of the parameters may be shared across the experiments can improve the identifiability of some compartmental models.

Skipping a pillar does not make for strong foundations: Pharmacokinetic-pharmacodynamic reasoning behind the shape of dose-response relationships in oncology.

Dose-response analysis is often applied to the quantification of drug-effect especially for slowly responding disease end points where a comparison is made across dose levels after a particular period of treatment. It has long been recognized that exposure - response is more appropriate than dose-response. However, trials necessarily are designed as dose-response experiments. Second, a wide range of functional forms are used to express relationships between dose and response. These considerations are also important for clinical development because pharmacokinetic (PK; and variability) plus pharmacokinetic-pharmacodynamic modeling may allow one to anticipate the shape of the dose-response curve and so the trial design. Here, we describe how the location and steepness of the dose response is determined by the PKs of the compound being tested and its exposure-response relationship in terms of potency (location), efficacy (maximum effect) and Hill coefficient (steepness). Thus, the location (50% effective dose [ED50 ]) is dependent not only on the potency (half-maximal effective concentration) but also the compound's PKs. Similarly, the steepness of the dose response is shown to be a function of the half-life of the drug. It is also shown that the shape of relationship varies dependent on the assumed time course of the disease. This is important in the context of drug-discovery where the in vivo potencies of compounds are compared as well as when considering an analysis of summary data (for example, model-based meta-analysis) for clinical decision making.

Population pharmacokinetic/pharmacodynamic modelling of the anti-TNF-α polyclonal fragment antibody AZD9773 in patients with severe sepsis.

AZD9773 is an ovine-derived, polyclonal, anti-tumour necrosis factor-alpha (TNF-α) antibody fragment. Using data from an AZD9773 Phase IIa study in patients with severe sepsis (clinicaltrials.gov: NCT00615017), a population pharmacokinetic/pharmacodynamic (PK/PD) model was developed. The model assessed the influence of various covariates on the PK of AZD9773 and the relationship between AZD9773 exposure and serological TNF-α concentration. A linear two-compartment model was used to describe AZD9773 concentration-time data. A stepwise covariate analysis was performed on the PK parameters. Subsequently, the serological TNF-α concentrations and drug effect were captured using an indirect response model, with a variable production rate of TNF-α described by a quadratic function. Creatinine clearance (CrCL) was the only covariate with a significant effect on the PK of AZD9773. A typical patient's drug clearance varied with CrCL; the relationship was non-linear. Diagnostic analysis of the PK/PD model showed that the fit was good, both across cohorts and in AZD9773-treated versus placebo patients. Serological TNF-α concentrations and the reduction of measurable serum TNF-α by AZD9773 were well characterized across all the cohorts evaluated in the Phase IIa study. This population PK/PD model was subsequently used to simulate alternative dosing options for a Phase IIb study (clinicaltrials.gov: NCT01145560).

Comparison of classical tumour growth models for patient derived and cell-line derived xenografts using the nonlinear mixed-effects framework.

In this study we compare seven mathematical models of tumour growth using nonlinear mixed-effects which allows for a simultaneous fitting of multiple data and an estimation of both mean behaviour and variability. This is performed for two large datasets, a patient-derived xenograft (PDX) dataset consisting of 220 PDXs spanning six different tumour types and a cell-line derived xenograft (CDX) dataset consisting of 25 cell lines spanning eight tumour types. Comparison of the models is performed by means of visual predictive checks (VPCs) as well as the Akaike Information Criterion (AIC). Additionally, we fit the models to 500 bootstrap samples drawn from the datasets to expand the comparison of the models under dataset perturbations and understand the growth kinetics that are best fitted by each model. Through qualitative and quantitative metrics the best models are identified the effectiveness and practicality of simpler models is highlighted.

The ingredients for an antimicrobial mathematical modelling broth.

Mathematical modelling has made significant contributions to the optimization of the use of antimicrobial treatments. This article discusses the key processes that such mathematical modelling should attempt to capture. In particular, this article highlights that the response of the host immune system requires quantification, and this is illustrated with a novel model structure.

How translational modeling in oncology needs to get the mechanism just right.

Translational model-based approaches have played a role in increasing success in the development of novel anticancer treatments. However, despite this, significant translational uncertainty remains from animal models to patients. Optimization of dose and scheduling (regimen) of drugs to maximize the therapeutic utility (maximize efficacy while avoiding limiting toxicities) is still predominately driven by clinical investigations. Here, we argue that utilizing pragmatic mechanism-based translational modeling of nonclinical data can further inform this optimization. Consequently, a prototype model is demonstrated that addresses the required fundamental mechanisms.

Step-by-step comparison of ordinary differential equation and agent-based approaches to pharmacokinetic-pharmacodynamic models.

Mathematical models in oncology aid in the design of drugs and understanding of their mechanisms of action by simulation of drug biodistribution, drug effects, and interaction between tumor and healthy cells. The traditional approach in pharmacometrics is to develop and validate ordinary differential equation models to quantify trends at the population level. In this approach, time-course of biological measurements is modeled continuously, assuming a homogenous population. Another approach, agent-based models, focuses on the behavior and fate of biological entities at the individual level, which subsequently could be summarized to reflect the population level. Heterogeneous cell populations and discrete events are simulated, and spatial distribution can be incorporated. In this tutorial, an agent-based model is presented and compared to an ordinary differential equation model for a tumor efficacy model inhibiting the pERK pathway. We highlight strengths, weaknesses, and opportunities of each approach.

Epstein-Barr Virus Reactivation After Paediatric Haematopoietic Stem Cell Transplantation: Risk Factors and Sensitivity Analysis of Mathematical Model.

Epstein-Barr virus (EBV) establishes a lifelong latent infection in healthy humans, kept under immune control by cytotoxic T cells (CTLs). Following paediatric haematopoetic stem cell transplantation (HSCT), a loss of immune surveillance leads to opportunistic outgrowth of EBV-infected cells, resulting in EBV reactivation, which can ultimately progress to post-transplant lymphoproliferative disorder (PTLD). The aims of this study were to identify risk factors for EBV reactivation in children in the first 100 days post-HSCT and to assess the suitability of a previously reported mathematical model to mechanistically model EBV reactivation kinetics in this cohort. Retrospective electronic data were collected from 56 children who underwent HSCT at Great Ormond Street Hospital (GOSH) between 2005 and 2016. Using EBV viral load (VL) measurements from weekly quantitative PCR (qPCR) monitoring post-HSCT, a multivariable Cox proportional hazards (Cox-PH) model was developed to assess time to first EBV reactivation event in the first 100 days post-HSCT. Sensitivity analysis of a previously reported mathematical model was performed to identify key parameters affecting EBV VL. Cox-PH modelling revealed EBV seropositivity of the HSCT recipient and administration of anti-thymocyte globulin (ATG) pre-HSCT to be significantly associated with an increased risk of EBV reactivation in the first 100 days post-HSCT (adjusted hazard ratio (AHR) = 2.32, P = 0.02; AHR = 2.55, P = 0.04). Five parameters were found to affect EBV VL in sensitivity analysis of the previously reported mathematical model. In conclusion, we have assessed the effect of multiple covariates on EBV reactivation in the first 100 days post-HSCT in children and have identified key parameters in a previously reported mechanistic mathematical model that affect EBV VL. Future work will aim to fit this model to patient EBV VLs, develop the model to account for interindividual variability and model the effect of clinically relevant covariates such as rituximab therapy and ATG on EBV VL.

A meta-analysis of tumour response and relapse kinetics based on 34,881 patients: A question of cancer type, treatment and line of treatment.

PURPOSE: Cancer disease burden is commonly assessed radiologically in solid tumours in support of response assessment via the RECIST criteria. These longitudinal data are amenable to mathematical modelling and these models characterise the initial tumour size, initial tumour shrinkage in responding patients and rate of regrowth as patient's disease progresses. Knowing how these parameters vary between patient populations and treatments would inform translational modelling approaches from non-clinical data as well as clinical trial design. EXPERIMENTAL DESIGN: Here a meta-analysis of reported model parameter values is reported. Appropriate literature was identified via a PubMed search and the application of text-based clustering approaches. The resulting parameter estimates are examined graphically and with ANOVA. RESULTS: Parameter values from a total of 80 treatment arms were identified based on 80 trial arms containing a total of 34,881 patients. Parameter estimates are generally consistent. It is found that a significant proportion of the variation in rates of tumour shrinkage and regrowth are explained by differing cancer and treatment: cancer type accounts for 66% of the variation in shrinkage rate and 71% of the variation in reported regrowth rates. Mean average parameter values by cancer and treatment are also reported. CONCLUSIONS: Mathematical modelling of longitudinal data is most often reported on a per clinical trial basis. However, the results reported here suggest that a more integrative approach would benefit the development of new treatments as well as the further optimisation of those currently used.

Preclinical Comparison of the Blood-brain barrier Permeability of Osimertinib with Other EGFR TKIs.

PURPOSE: Osimertinib is a potent and selective EGFR tyrosine kinase inhibitor (EGFR-TKI) of both sensitizing and T790M resistance mutations. To treat metastatic brain disease, blood-brain barrier (BBB) permeability is considered desirable for increasing clinical efficacy. EXPERIMENTAL DESIGN: We examined the level of brain penetration for 16 irreversible and reversible EGFR-TKIs using multiple in vitro and in vivo BBB preclinical models. RESULTS: In vitro osimertinib was the weakest substrate for human BBB efflux transporters (efflux ratio 3.2). In vivo rat free brain to free plasma ratios (Kpuu) show osimertinib has the most BBB penetrance (0.21), compared with the other TKIs (Kpuu ≤ 0.12). PET imaging in Cynomolgus macaques demonstrated osimertinib was the only TKI among those tested to achieve significant brain penetrance (C max %ID 1.5, brain/blood Kp 2.6). Desorption electrospray ionization mass spectroscopy images of brains from mouse PC9 macrometastases models showed osimertinib readily distributes across both healthy brain and tumor tissue. Comparison of osimertinib with the poorly BBB penetrant afatinib in a mouse PC9 model of subclinical brain metastases showed only osimertinib has a significant effect on rate of brain tumor growth. CONCLUSIONS: These preclinical studies indicate that osimertinib can achieve significant exposure in the brain compared with the other EGFR-TKIs tested and supports the ongoing clinical evaluation of osimertinib for the treatment of EGFR-mutant brain metastasis. This work also demonstrates the link between low in vitro transporter efflux ratios and increased brain penetrance in vivo supporting the use of in vitro transporter assays as an early screen in drug discovery.

Clone Wars: Quantitatively Understanding Cancer Drug Resistance.

A key aim of early clinical development for new cancer treatments is to detect the potential for efficacy early and to identify a safe therapeutic dose to take forward to phase II. Because of this need, researchers have sought to build mathematical models linking initial radiologic tumor response, often assessed after 6 to 8 weeks of treatment, with overall survival. However, there has been mixed success of this approach in the literature. We argue that evolutionary selection pressure should be considered to interpret these early efficacy signals and so optimize cancer therapy.