Advancing cancer drug development with mechanistic mathematical modeling: bridging the gap between theory and practice.

Quantitative predictive modeling of cancer growth, progression, and individual response to therapy is a rapidly growing field. Researchers from mathematical modeling, systems biology, pharmaceutical industry, and regulatory bodies, are collaboratively working on predictive models that could be applied for drug development and, ultimately, the clinical management of cancer patients. A plethora of modeling paradigms and approaches have emerged, making it challenging to compile a comprehensive review across all subdisciplines. It is therefore critical to gauge fundamental design aspects against requirements, and weigh opportunities and limitations of the different model types. In this review, we discuss three fundamental types of cancer models: space-structured models, ecological models, and immune system focused models. For each type, it is our goal to illustrate which mechanisms contribute to variability and heterogeneity in cancer growth and response, so that the appropriate architecture and complexity of a new model becomes clearer. We present the main features addressed by each of the three exemplary modeling types through a subjective collection of literature and illustrative exercises to facilitate inspiration and exchange, with a focus on providing a didactic rather than exhaustive overview. We close by imagining a future multi-scale model design to impact critical decisions in oncology drug development.

An asymptotic description of a basic FcRn-regulated clearance mechanism and its implications for PBPK modelling of large antibodies.

A basic FcRn-regulated clearance mechanism is investigated using the method of matched asymptotic expansions. The broader aim of the work is to obtain further insight on the mechanism, thereby providing theoretical support for future pharmacologically-based pharmacokinetic modelling efforts. The corresponding governing equations are first non-dimensionalised and the order of magnitudes of the model parameters are assessed based on their values reported in the literature. Under the assumption of high FcRn-binding affinity, analytical approximations are derived that are valid over the characteristic phases of the problem. Additionally, relatively simple equations relating clearance and AUC to physiological model parameters are derived, which are valid over the longest characteristic time scale of the problem. For lower to moderate doses clearance is effectively linear, whereas for higher doses it is nonlinear. It is shown that for all doses sufficiently high the leading-order approximation for the IgG concentration in plasma, over the longest characteristic time scale, is independent of the initial dose. This is because IgG that is in 'excess' of FcRn is eliminated over a time scale much shorter than that of the terminal phase. In conclusion, analytical approximations of the basic FcRn mechanism have been derived using matched asymptotic expansions, leading to a simple equation relating clearance to FcRn binding affinity, the ratio of degradation and FcRn concentration, and the volumes of the system.

Assessing the performance of QSP models: biology as the driver for validation.

Validation of a quantitative model is a critical step in establishing confidence in the model's suitability for whatever analysis it was designed. While processes for validation are well-established in the statistical sciences, the field of quantitative systems pharmacology (QSP) has taken a more piecemeal approach to defining and demonstrating validation. Although classical statistical methods can be used in a QSP context, proper validation of a mechanistic systems model requires a more nuanced approach to what precisely is being validated, and what role said validation plays in the larger context of the analysis. In this review, we summarize current thoughts of QSP validation in the scientific community, contrast the aims of statistical validation from several contexts (including inference, pharmacometrics analysis, and machine learning) with the challenges faced in QSP analysis, and use examples from published QSP models to define different stages or levels of validation, any of which may be sufficient depending on the context at hand.

Evaluating a Multiscale Mechanistic Model of the Immune System to Predict Human Immunogenicity for a Biotherapeutic in Phase 1.

A mechanistic model of the immune response was evaluated for its ability to predict anti-drug antibody (ADA) and their impact on pharmacokinetics (PK) and pharmacodynamics (PD) for a biotherapeutic in a phase 1 clinical trial. Observed ADA incidence ranged from 33 to 67% after single doses and 27-50% after multiple doses. The model captured the single dose incidence well; however, there was overprediction after multiple dosing. The model was updated to include a T-regulatory (Treg) cell mediated tolerance, which reduced the overprediction (relative decrease in predicted incidence rate of 21.5-59.3% across multidose panels) without compromising the single dose predictions (relative decrease in predicted incidence rate of 0.6-13%). The Treg-adjusted model predicted no ADA impact on PK or PD, consistent with the observed data. A prospective phase 2 trial was simulated, including co-medication effects in the form of corticosteroid-induced immunosuppression. Predicted ADA incidences were 0-10%, depending on co-medication dosage. This work demonstrates the utility in applying an integrated, iterative modeling approach to predict ADA during different stages of clinical development.

A Survey of Software Tool Utilization and Capabilities for Quantitative Systems Pharmacology: What We Have and What We Need.

Quantitative systems pharmacology (QSP) is a rapidly emerging discipline with application across a spectrum of challenges facing the pharmaceutical industry, including mechanistically informed prioritization of target pathways and combinations in discovery, target population, and dose expansion decisions early in clinical development, and analyses for regulatory authorities late in clinical development. QSP's development has influences from physiologic modeling, systems biology, physiologically-based pharmacokinetic modeling, and pharmacometrics. Given a varied scientific heritage, a variety of tools to accomplish the demands of model development, application, and model-based analysis of available data have been developed. We report the outcome from a community survey and resulting analysis of how modelers view the impact and growth of QSP, how they utilize existing tools, and capabilities they need improved to further accelerate their impact on drug development. These results serve as a benchmark and roadmap for advancements to the QSP tool set.

QSP Toolbox: Computational Implementation of Integrated Workflow Components for Deploying Multi-Scale Mechanistic Models.

Quantitative systems pharmacology (QSP) modeling has become increasingly important in pharmaceutical research and development, and is a powerful tool to gain mechanistic insights into the complex dynamics of biological systems in response to drug treatment. However, even once a suitable mathematical framework to describe the pathophysiology and mechanisms of interest is established, final model calibration and the exploration of variability can be challenging and time consuming. QSP models are often formulated as multi-scale, multi-compartment nonlinear systems of ordinary differential equations. Commonly accepted modeling strategies, workflows, and tools have promise to greatly improve the efficiency of QSP methods and improve productivity. In this paper, we present the QSP Toolbox, a set of functions, structure array conventions, and class definitions that computationally implement critical elements of QSP workflows including data integration, model calibration, and variability exploration. We present the application of the toolbox to an ordinary differential equations-based model for antibody drug conjugates. As opposed to a single stepwise reference model calibration, the toolbox also facilitates simultaneous parameter optimization and variation across multiple in vitro, in vivo, and clinical assays to more comprehensively generate alternate mechanistic hypotheses that are in quantitative agreement with available data. The toolbox also includes scripts for developing and applying virtual populations to mechanistic exploration of biomarkers and efficacy. We anticipate that the QSP Toolbox will be a useful resource that will facilitate implementation, evaluation, and sharing of new methodologies in a common framework that will greatly benefit the community.

Alzheimer's disease: rapid and slow progression.

The variability in the progression of Alzheimer's disease (AD) across patients has made identification of disease-delaying treatments difficult. Quantitative analysis of this variability has important implications in understanding the pathophysiology of AD and identifying disease-delaying treatments. The functional assessment staging (FAST) procedure characterizes seven stages in the course of AD from normal ageing to severe dementia. The present study applied statistical methods to analyse FAST stage durations from a dataset of 648 AD patients. These methods uncovered two distinct types of disease progression, characterized by different mean progression rates. We identified two separate distributions of FAST stage progression times differing by up to 2 years in mean duration within each stage. These results further indicate that if a patient progresses rapidly through a given FAST stage, then their further progression is also likely to be rapid. These findings support the hypothesis that progression of AD can occur via two different pathophysiological mechanisms that lead to distinct average rates of decline.

Calculating stage duration statistics in multistage diseases.

Many human diseases are characterized by multiple stages of progression. While the typical sequence of disease progression can be identified, there may be large individual variations among patients. Identifying mean stage durations and their variations is critical for statistical hypothesis testing needed to determine if treatment is having a significant effect on the progression, or if a new therapy is showing a delay of progression through a multistage disease. In this paper we focus on two methods for extracting stage duration statistics from longitudinal datasets: an extension of the linear regression technique, and a counting algorithm. Both are non-iterative, non-parametric and computationally cheap methods, which makes them invaluable tools for studying the epidemiology of diseases, with a goal of identifying different patterns of progression by using bioinformatics methodologies. Here we show that the regression method performs well for calculating the mean stage durations under a wide variety of assumptions, however, its generalization to variance calculations fails under realistic assumptions about the data collection procedure. On the other hand, the counting method yields reliable estimations for both means and variances of stage durations. Applications to Alzheimer disease progression are discussed.

High degree of heterogeneity in Alzheimer's disease progression patterns.

There have been several reports on the varying rates of progression among Alzheimer's Disease (AD) patients; however, there has been no quantitative study of the amount of heterogeneity in AD. Obtaining a reliable quantitative measure of AD progression rates and their variances among the patients for each stage of AD is essential for evaluating results of any clinical study. The Global Deterioration Scale (GDS) and Functional Assessment Staging procedure (FAST) characterize seven stages in the course of AD from normal aging to severe dementia. Each GDS/FAST stage has a published mean duration, but the variance is unknown. We use statistical analysis to reconstruct GDS/FAST stage durations in a cohort of 648 AD patients with an average follow-up time of 4.78 years. Calculations for GDS/FAST stages 4-6 reveal that the standard deviations for stage durations are comparable with their mean values, indicating the presence of large variations in the AD progression among patients. Such amount of heterogeneity in the course of progression of AD is consistent with the existence of several sub-groups of AD patients, which differ by their patterns of decline.

Signal response sensitivity in the yeast mitogen-activated protein kinase cascade.

The yeast pheromone response pathway is a canonical three-step mitogen activated protein kinase (MAPK) cascade which requires a scaffold protein for proper signal transduction. Recent experimental studies into the role the scaffold plays in modulating the character of the transduced signal, show that the presence of the scaffold increases the biphasic nature of the signal response. This runs contrary to prior theoretical investigations into how scaffolds function. We describe a mathematical model of the yeast MAPK cascade specifically designed to capture the experimental conditions and results of these empirical studies. We demonstrate how the system can exhibit either graded or ultrasensitive (biphasic) response dynamics based on the binding kinetics of enzymes to the scaffold. At the basis of our theory is an analytical result that weak interactions make the response biphasic while tight interactions lead to a graded response. We then show via an analysis of the kinetic binding rate constants how the results of experimental manipulations, modeled as changes to certain of these binding constants, lead to predictions of pathway output consistent with experimental observations. We demonstrate how the results of these experimental manipulations are consistent within the framework of our theoretical treatment of this scaffold-dependent MAPK cascades, and how future efforts in this style of systems biology can be used to interpret the results of other signal transduction observations.

Selection in spatial stochastic models of cancer: migration as a key modulator of fitness.

BACKGROUND: We study the selection dynamics in a heterogeneous spatial colony of cells. We use two spatial generalizations of the Moran process, which include cell divisions, death and migration. In the first model, migration is included explicitly as movement to a proximal location. In the second, migration is implicit, through the varied ability of cell types to place their offspring a distance away, in response to another cell's death. RESULTS: In both models, we find that migration has a direct positive impact on the ability of a single mutant cell to invade a pre-existing colony. Thus, a decrease in the growth potential can be compensated by an increase in cell migration. We further find that the neutral ridges (the set of all types with the invasion probability equal to that of the host cells) remain invariant under the increase of system size (for large system sizes), thus making the invasion probability a universal characteristic of the cells selection status. We find that repeated instances of large scale cell-death, such as might arise during therapeutic intervention or host response, strongly select for the migratory phenotype. CONCLUSIONS: These models can help explain the many examples in the biological literature, where genes involved in cell's migratory and invasive machinery are also associated with increased cellular fitness, even though there is no known direct effect of these genes on the cellular reproduction. The models can also help to explain how chemotherapy may provide a selection mechanism for highly invasive phenotypes.

Specificity and robustness of the mammalian MAPK-IEG network.

The mitogen-activated protein kinase cascade is a conserved signal transduction pathway found in organisms of complexity spanning from yeast to humans. In many mammalian tissue types, this pathway can correctly transduce signals from different extracellular messengers, leading to specific and often mutually exclusive cellular responses. The transduced signal is tuned by a complicated set of positive and negative feedback control mechanisms and fed into a downstream gene expression network. This network, based on the immediate early gene system, has two possible, mutually exclusive outcomes. Using a mathematical model, we study how different stimuli lead to different temporal signal structure. Further, we investigate how each of the feedback controls contributes to the overall specificity of the gene expression output, and hypothesize that the complicated nature of the mammalian mitogen-activated protein kinase pathway results in a system able to robustly identify and transduce the proper signal without investing in two completely separate signal cascades. Finally, we quantify the role of the RKIP protein in shaping the signal, and propose a novel mechanism of its involvement in cancer metastasis.

Explicit separation of growth and motility in a new tumor cord model.

We investigate a new model of tumor growth in which cell motility is considered an explicitly separate process from growth. Bulk tumor expansion is modeled by individual cell motility in a density-dependent diffusion process. This model is implemented in the context of an in vivo system, the tumor cord. We investigate numerically microscale density distributions of different cell classes and macroscale whole tumor growth rates as functions of the strength of transitions between classes. Our results indicate that the total tumor growth follows a classical von Bertalanffy growth profile, as many in vivo tumors are observed to do. This provides a quick validation for the model hypotheses. The microscale and macroscale properties are both sensitive to fluctuations in the transition parameters, and grossly adopt one of two phenotypic profiles based on their parameter regime. We analyze these profiles and use the observations to classify parameter regimes by their phenotypes. This classification yields a novel hypothesis for the early evolutionary selection of the metastatic phenotype by selecting against less motile cells which grow to higher densities and may therefore induce local collapse of the vascular network.

Identification of histidine tautomers in proteins by 2D 1H/13C(delta2) one-bond correlated NMR.

If the 13Cdelta2 chemical shift of neutral ("high pH") histidine is >122 ppm, primarily Ndelta1-H tautomer (2) is indicated; if it is <122 ppm, primarily Nepsilon2-H tautomer (1) is indicated. His resonances from the catalytic triad of active serine proteases, for example, are readily distinguished from those of denatured enzyme. The 13Cdelta2 chemical shifts increased by 6.2 ppm for the catalytic histidines in both alpha-lytic protease and subtilisin BPN' in raising the pH from that of imidazolium cation to that of tautomer 2. This tautomer identification method is easy to implement, requiring only bioincorporation of [U-13C] (or the more readily available [U-13C,15N])-histidine. Standard 1H/13C correlation HMQC or HSQC NMR pulse programs then yield the 13Cdelta2 chemical shifts with the benefit of high 1H sensitivity. Because of large one-bond spin-couplings (1JCH approximately 200 Hz), the method should extend to proteins having large 1H and 13C line widths, including very high molecular weights.

Tautomerism, acid-base equilibria, and H-bonding of the six histidines in subtilisin BPN' by NMR.

We have determined by (15)N, (1)H, and (13)C NMR, the chemical behavior of the six histidines in subtilisin BPN' and their PMSF and peptide boronic acid complexes in aqueous solution as a function of pH in the range of from 5 to 11, and have assigned every (15)N, (1)H, C(epsilon 1), and C(delta2) resonance of all His side chains in resting enzyme. Four of the six histidine residues (17, 39, 67, and 226) are neutrally charged and do not titrate. One histidine (238), located on the protein surface, titrates with pK(a) = 7.30 +/- 0.03 at 25 degrees C, having rapid proton exchange, but restricted mobility. The active site histidine (64) in mutant N155A titrates with a pK(a) value of 7.9 +/- 0.3 and sluggish proton exchange behavior, as shown by two-site exchange computer lineshape simulation. His 64 in resting enzyme contains an extremely high C(epsilon 1)-H proton chemical shift of 9.30 parts per million (ppm) owing to a conserved C(epsilon 1)-H(.)O=C H-bond from the active site imidazole to a backbone carbonyl group, which is found in all known serine proteases representing all four superfamilies. Only His 226, and His 64 at high pH, exist as the rare N(delta1)-H tautomer, exhibiting (13)C(delta1) chemical shifts approximately 9 ppm higher than those for N(epsilon 2)-H tautomers. His 64 in the PMSF complex, unlike that in the resting enzyme, is highly mobile in its low pH form, as shown by (15)N-(1)H NOE effects, and titrates with rapid proton exchange kinetics linked to a pK(a) value of 7.47 +/- 0.02.