Connecting Agent-Based Models with High-Dimensional Parameter Spaces to Multidimensional Data Using SMoRe ParS: A Surrogate Modeling Approach.

Across a broad range of disciplines, agent-based models (ABMs) are increasingly utilized for replicating, predicting, and understanding complex systems and their emergent behavior. In the biological and biomedical sciences, researchers employ ABMs to elucidate complex cellular and molecular interactions across multiple scales under varying conditions. Data generated at these multiple scales, however, presents a computational challenge for robust analysis with ABMs. Indeed, calibrating ABMs remains an open topic of research due to their own high-dimensional parameter spaces. In response to these challenges, we extend and validate our novel methodology, Surrogate Modeling for Reconstructing Parameter Surfaces (SMoRe ParS), arriving at a computationally efficient framework for connecting high dimensional ABM parameter spaces with multidimensional data. Specifically, we modify SMoRe ParS to initially confine high dimensional ABM parameter spaces using unidimensional data, namely, single time-course information of in vitro cancer cell growth assays. Subsequently, we broaden the scope of our approach to encompass more complex ABMs and constrain parameter spaces using multidimensional data. We explore this extension with in vitro cancer cell inhibition assays involving the chemotherapeutic agent oxaliplatin. For each scenario, we validate and evaluate the effectiveness of our approach by comparing how well ABM simulations match the experimental data when using SMoRe ParS-inferred parameters versus parameters inferred by a commonly used direct method. In so doing, we show that our approach of using an explicitly formulated surrogate model as an interlocutor between the ABM and the experimental data effectively calibrates the ABM parameter space to multidimensional data. Our method thus provides a robust and scalable strategy for leveraging multidimensional data to inform multiscale ABMs and explore the uncertainty in their parameters.

Understanding the mechanisms of HPV-related carcinogenesis: Implications for cell cycle dynamics.

The role of human papillomavirus (HPV) as a causative agent for epithelial cancers is well-known, but many open questions remain regarding the downstream gene regulatory effects of viral proteins E6 and E7 on the cell cycle. Here, we extend a cell cycle model originally presented by Gérard and Goldbeter (2009) in order to capture the effects of E6 and E7 on key actors in the cell cycle. Results suggest that E6 is sufficient to reverse p53-induced quiescence, while E7 is sufficient to reverse p16INK4a-induced quiescence; both E6 and E7 are necessary when p53 and p16INK4a are both active. Moreover, E7 appears to play a role as a "growth factor substitute", inducing cell division in the absence of growth factor. Low levels of E7 may permit regular cell division, but the results suggest that higher levels of E7 dysregulate the cell cycle in ways that may destabilize the cellular genome. The mechanisms explored here provide opportunities for developing new treatment targets that take advantage of the cell cycle regulatory system to prevent HPV-related cancer effects.

Tumor growth towards lower extracellular matrix conductivity regions under Darcy's Law and steady morphology.

We study a classic Darcy's law model for tumor cell motion with inhomogeneous and isotropic conductivity. The tumor cells are assumed to be a constant density fluid flowing through porous extracellular matrix (ECM). The ECM is assumed to be rigid and motionless with constant porosity. One and two dimensional simulations show that the tumor mass grows from high to low conductivity regions when the tumor morphology is steady. In the one-dimensional case, we proved that when the tumor size is steady, the tumor grows towards lower conductivity regions. We conclude that this phenomenon is produced by the coupling of a special inward flow pattern in the steady tumor and Darcy's law which gives faster flow speed in higher conductivity regions.

A mathematical model for IL-6-mediated, stem cell driven tumor growth and targeted treatment.

Targeting key regulators of the cancer stem cell phenotype to overcome their critical influence on tumor growth is a promising new strategy for cancer treatment. Here we present a modeling framework that operates at both the cellular and molecular levels, for investigating IL-6 mediated, cancer stem cell driven tumor growth and targeted treatment with anti-IL6 antibodies. Our immediate goal is to quantify the influence of IL-6 on cancer stem cell self-renewal and survival, and to characterize the subsequent impact on tumor growth dynamics. By including the molecular details of IL-6 binding, we are able to quantify the temporal changes in fractional occupancies of bound receptors and their influence on tumor volume. There is a strong correlation between the model output and experimental data for primary tumor xenografts. We also used the model to predict tumor response to administration of the humanized IL-6R monoclonal antibody, tocilizumab (TCZ), and we found that as little as 1mg/kg of TCZ administered weekly for 7 weeks is sufficient to result in tumor reduction and a sustained deceleration of tumor growth.

Mathematical Modeling of Cellular Cross-Talk Between Endothelial and Tumor Cells Highlights Counterintuitive Effects of VEGF-Targeted Therapies.

Tumor growth and progression are critically dependent on the establishment of a vascular support system. This is often accomplished via the expression of pro-angiogenic growth factors, including members of the vascular endothelial growth factor (VEGF) family of ligands. VEGF ligands are overexpressed in a wide variety of solid tumors and therefore have inspired optimism that inhibition of the different axes of the VEGF pathway-alone or in combination-would represent powerful anti-angiogenic therapies for most cancer types. When considering treatments that target VEGF and its receptors, it is difficult to tease out the differential anti-angiogenic and anti-tumor effects of all combinations experimentally because tumor cells and vascular endothelial cells are engaged in a dynamic cross-talk that impacts key aspects of tumorigenesis, independent of angiogenesis. Here we develop a mathematical model that connects intracellular signaling responsible for both endothelial and tumor cell proliferation and death to population-level cancer growth and angiogenesis. We use this model to investigate the effect of bidirectional communication between endothelial cells and tumor cells on treatments targeting VEGF and its receptors both in vitro and in vivo. Our results underscore the fact that in vitro therapeutic outcomes do not always translate to the in vivo situation. For example, our model predicts that certain therapeutic combinations result in antagonism in vivo that is not observed in vitro. Mathematical modeling in this direction can shed light on the mechanisms behind experimental observations that manipulating VEGF and its receptors is successful in some cases but disappointing in others.

Modeling head and neck cancer stem cell-mediated tumorigenesis.

A large body of literature has emerged supporting the importance of cancer stem cells (CSCs) in the pathogenesis of head and neck cancers. CSCs are a subpopulation of cells within a tumor that share the properties of self-renewal and multipotency with stem cells from normal tissue. Their functional relevance to the pathobiology of cancer arises from the unique properties of tumorigenicity, chemotherapy resistance, and their ability to metastasize and invade distant tissues. Several molecular profiles have been used to discriminate a stem cell from a non-stem cell. CSCs can be grown for study and further enriched using a number of in vitro techniques. An evolving option for translational research is the use of mathematical and computational models to describe the role of CSCs in complex tumor environments. This review is focused discussing the evidence emerging from modeling approaches that have clarified the impact of CSCs to the biology of cancer.

A phase II trial of the BCL-2 homolog domain 3 mimetic AT-101 in combination with docetaxel for recurrent, locally advanced, or metastatic head and neck cancer.

BACKGROUND: AT-101 is a BCL-2 Homolog domain 3 mimetic previously demonstrated to have tumoricidal effects in advanced solid organ malignancies. Given the evidence of activity in xenograft models, treatment with AT-101 in combination with docetaxel is a therapeutic doublet of interest in metastatic head and neck squamous cell carcinoma. PATIENTS AND METHODS: Patients included in this trial had unresectable, recurrent, or distantly metastatic head and neck squamous cell carcinoma (R/M HNSCC) not amenable to curative radiation or surgery. This was an open label randomized, phase II trial in which patients were administered AT-101 in addition to docetaxel. The three treatment arms were docetaxel, docetaxel plus pulse dose AT-101, and docetaxel plus metronomic dose AT-101. The primary endpoint of this trial was overall response rate. RESULTS: Thirty-five patients were registered and 32 were evaluable for treatment response. Doublet therapy with AT-101 and docetaxel was well tolerated with only 2 patients discontinuing therapy due to treatment related toxicities. The overall response rate was 11 % (4 partial responses) with a clinical benefit rate of 74 %. Median progression free survival was 4.3 months (range: 0.7-13.7) and overall survival was 5.5 months (range: 0.4-24). No significant differences were noted between dosing strategies. CONCLUSION: Although met with a favorable toxicity profile, the addition of AT-101 to docetaxel in R/M HNSCC does not appear to demonstrate evidence of efficacy.

Dysregulated FGFR3 signaling alters the immune landscape in bladder cancer and presents therapeutic possibilities in an agent-based model.

Bladder cancer is an increasingly prevalent global disease that continues to cause morbidity and mortality despite recent advances in treatment. Immune checkpoint inhibitors (ICI) and fibroblast growth factor receptor (FGFR)-targeted therapeutics have had modest success in bladder cancer when used as monotherapy. Emerging data suggests that the combination of these two therapies could lead to improved clinical outcomes, but the optimal strategy for combining these agents remains uncertain. Mathematical models, specifically agent-based models (ABMs), have shown recent successes in uncovering the multiscale dynamics that shape the trajectory of cancer. They have enabled the optimization of treatment methods and the identification of novel therapeutic strategies. To assess the combined effects of anti-PD-1 and anti-FGFR3 small molecule inhibitors (SMI) on tumor growth and the immune response, we built an ABM that captures key facets of tumor heterogeneity and CD8+ T cell phenotypes, their spatial interactions, and their response to therapeutic pressures. Our model quantifies how tumor antigenicity and FGFR3 activating mutations impact disease trajectory and response to anti-PD-1 antibodies and anti-FGFR3 SMI. We find that even a small population of weakly antigenic tumor cells bearing an FGFR3 mutation can render the tumor resistant to combination therapy. However, highly antigenic tumors can overcome therapeutic resistance mediated by FGFR3 mutation. The optimal therapy depends on the strength of the FGFR3 signaling pathway. Under certain conditions, ICI alone is optimal; in others, ICI followed by anti-FGFR3 therapy is best. These results indicate the need to quantify FGFR3 signaling and the fitness advantage conferred on bladder cancer cells harboring this mutation. This ABM approach may enable rationally designed treatment plans to improve clinical outcomes.

Haematopoietic stem cells do not asymmetrically segregate chromosomes or retain BrdU.

Stem cells are proposed to segregate chromosomes asymmetrically during self-renewing divisions so that older ('immortal') DNA strands are retained in daughter stem cells whereas newly synthesized strands segregate to differentiating cells. Stem cells are also proposed to retain DNA labels, such as 5-bromo-2-deoxyuridine (BrdU), either because they segregate chromosomes asymmetrically or because they divide slowly. However, the purity of stem cells among BrdU-label-retaining cells has not been documented in any tissue, and the 'immortal strand hypothesis' has not been tested in a system with definitive stem cell markers. Here we tested these hypotheses in haematopoietic stem cells (HSCs), which can be highly purified using well characterized markers. We administered BrdU to newborn mice, mice treated with cyclophosphamide and granulocyte colony-stimulating factor, and normal adult mice for 4 to 10 days, followed by 70 days without BrdU. In each case, less than 6% of HSCs retained BrdU and less than 0.5% of all BrdU-retaining haematopoietic cells were HSCs, revealing that BrdU has poor specificity and poor sensitivity as an HSC marker. Sequential administration of 5-chloro-2-deoxyuridine and 5-iodo-2-deoxyuridine indicated that all HSCs segregate their chromosomes randomly. Division of individual HSCs in culture revealed no asymmetric segregation of the label. Thus, HSCs cannot be identified on the basis of BrdU-label retention and do not retain older DNA strands during division, indicating that these are not general properties of stem cells.

Phenotype switching in a global method for agent-based models of biological tissue.

Agent-based models (ABMs) are an increasingly important tool for understanding the complexities presented by phenotypic and spatial heterogeneity in biological tissue. The resolution a modeler can achieve in these regards is unrivaled by other approaches. However, this comes at a steep computational cost limiting either the scale of such models or the ability to explore, parameterize, analyze, and apply them. When the models involve molecular-level dynamics, especially cell-specific dynamics, the limitations are compounded. We have developed a global method for solving these computationally expensive dynamics significantly decreases the computational time without altering the behavior of the system. Here, we extend this method to the case where cells can switch phenotypes in response to signals in the microenvironment. We find that the global method in this context preserves the temporal population dynamics and the spatial arrangements of the cells while requiring markedly less simulation time. We thus add a tool for efficiently simulating ABMs that captures key facets of the molecular and cellular dynamics in heterogeneous tissue.

Mathematical model predicts tumor control patterns induced by fast and slow cytotoxic T lymphocyte killing mechanisms.

Immunotherapy has dramatically transformed the cancer treatment landscape largely due to the efficacy of immune checkpoint inhibitors (ICIs). Although ICIs have shown promising results for many patients, the low response rates in many cancers highlight the ongoing challenges in cancer treatment. Cytotoxic T lymphocytes (CTLs) execute their cell-killing function via two distinct mechanisms: a fast-acting, perforin-mediated process and a slower, Fas ligand (FasL)-driven pathway. Evidence also suggests that the preferred killing mechanism of CTLs depends on the antigenicity of tumor cells. To determine the critical factors affecting responses to ICIs, we construct an ordinary differential equation model describing in vivo tumor-immune dynamics in the presence of active or blocked PD-1/PD-L1 immune checkpoint. Specifically, we identify important aspects of the tumor-immune landscape that affect tumor size and composition in the short and long term. We also generate a virtual cohort of mice with diverse tumor and immune attributes to simulate the outcomes of immune checkpoint blockade in a heterogeneous population. By identifying key tumor and immune characteristics associated with tumor elimination, dormancy, and escape, we predict which fraction of a population potentially responds well to ICIs and ways to enhance therapeutic outcomes with combination therapy.

Implementing a global approach for efficiently simulating molecular dynamics in agent-based models of biological tissue.

This protocol explains how to take an agent-based model (ABM) with molecular dynamics and set it up to solve the molecular dynamics with a global approach. It can be used to speed up simulations significantly while retaining high levels of accuracy with the original ABM. Two options are presented for implementing this global approach, depending on the desired spatial variability in molecular concentrations. Both options coarse-grain the molecular dynamics in space by dividing the microenvironment into regions with uniform concentrations. For complete details on the use and execution of this protocol, please refer to Bergman et al. (2022).

Analysis of HIV latent infection model with multiple infection stages and different drug classes.

Latently infected CD4+ T cells represent one of the major obstacles to HIV eradication even after receiving prolonged highly active anti-retroviral therapy (HAART). Long-term use of HAART causes the emergence of drug-resistant virus which is then involved in HIV transmission. In this paper, we develop mathematical HIV models with staged disease progression by incorporating entry inhibitor and latently infected cells. We find that entry inhibitor has the same effect as protease inhibitor on the model dynamics and therefore would benefit HIV patients who developed resistance to many of current anti-HIV medications. Numerical simulations illustrate the theoretical results and show that the virus and latently infected cells reach an infected steady state in the absence of treatment and are eliminated under treatment whereas the model including homeostatic proliferation of latently infected cells maintains the virus at low level during suppressive treatment. Therefore, complete cure of HIV needs complete eradication of latent reservoirs.

A global method for fast simulations of molecular dynamics in multiscale agent-based models of biological tissues.

Agent-based models (ABMs) are a natural platform for capturing the multiple time and spatial scales in biological processes. However, these models are computationally expensive, especially when including molecular-level effects. The traditional approach to simulating this type of multiscale ABM is to solve a system of ordinary differential equations for the molecular events per cell. This significantly adds to the computational cost of simulations as the number of agents grows, which contributes to many ABMs being limited to around 10 5 cells. We propose an approach that requires the same computational time independent of the number of agents. This speeds up the entire simulation by orders of magnitude, allowing for more thorough explorations of ABMs with even larger numbers of agents. We use two systems to show that the new method strongly agrees with the traditionally used approach. This computational strategy can be applied to a wide range of biological investigations.

SMoRe ParS: A novel methodology for bridging modeling modalities and experimental data applied to 3D vascular tumor growth.

Multiscale systems biology is having an increasingly powerful impact on our understanding of the interconnected molecular, cellular, and microenvironmental drivers of tumor growth and the effects of novel drugs and drug combinations for cancer therapy. Agent-based models (ABMs) that treat cells as autonomous decision-makers, each with their own intrinsic characteristics, are a natural platform for capturing intratumoral heterogeneity. Agent-based models are also useful for integrating the multiple time and spatial scales associated with vascular tumor growth and response to treatment. Despite all their benefits, the computational costs of solving agent-based models escalate and become prohibitive when simulating millions of cells, making parameter exploration and model parameterization from experimental data very challenging. Moreover, such data are typically limited, coarse-grained and may lack any spatial resolution, compounding these challenges. We address these issues by developing a first-of-its-kind method that leverages explicitly formulated surrogate models (SMs) to bridge the current computational divide between agent-based models and experimental data. In our approach, Surrogate Modeling for Reconstructing Parameter Surfaces (SMoRe ParS), we quantify the uncertainty in the relationship between agent-based model inputs and surrogate model parameters, and between surrogate model parameters and experimental data. In this way, surrogate model parameters serve as intermediaries between agent-based model input and data, making it possible to use them for calibration and uncertainty quantification of agent-based model parameters that map directly onto an experimental data set. We illustrate the functionality and novelty of Surrogate Modeling for Reconstructing Parameter Surfaces by applying it to an agent-based model of 3D vascular tumor growth, and experimental data in the form of tumor volume time-courses. Our method is broadly applicable to situations where preserving underlying mechanistic information is of interest, and where computational complexity and sparse, noisy calibration data hinder model parameterization.

An Agent-Based Model of Combination Oncolytic Viral Therapy and Anti-PD-1 Immunotherapy Reveals the Importance of Spatial Location When Treating Glioblastoma.

Oncolytic viral therapies and immunotherapies are of growing clinical interest due to their selectivity for tumor cells over healthy cells and their immunostimulatory properties. These treatment modalities provide promising alternatives to the standard of care, particularly for cancers with poor prognoses, such as the lethal brain tumor glioblastoma (GBM). However, uncertainty remains regarding optimal dosing strategies, including how the spatial location of viral doses impacts therapeutic efficacy and tumor landscape characteristics that are most conducive to producing an effective immune response. We develop a three-dimensional agent-based model (ABM) of GBM undergoing treatment with a combination of an oncolytic Herpes Simplex Virus and an anti-PD-1 immunotherapy. We use a mechanistic approach to model the interactions between distinct populations of immune cells, incorporating both innate and adaptive immune responses to oncolytic viral therapy and including a mechanism of adaptive immune suppression via the PD-1/PD-L1 checkpoint pathway. We utilize the spatially explicit nature of the ABM to determine optimal viral dosing in both the temporal and spatial contexts. After proposing an adaptive viral dosing strategy that chooses to dose sites at the location of highest tumor cell density, we find that, in most cases, this adaptive strategy produces a more effective treatment outcome than repeatedly dosing in the center of the tumor.

Standing Variations Modeling Captures Inter-Individual Heterogeneity in a Deterministic Model of Prostate Cancer Response to Combination Therapy.

Sipuleucel-T (Provenge) is the first live cell vaccine approved for advanced, hormonally refractive prostate cancer. However, survival benefit is modest and the optimal combination or schedule of sipuleucel-T with androgen depletion remains unknown. We employ a nonlinear dynamical systems approach to modeling the response of hormonally refractive prostate cancer to sipuleucel-T. Our mechanistic model incorporates the immune response to the cancer elicited by vaccination, and the effect of androgen depletion therapy. Because only a fraction of patients benefit from sipuleucel-T treatment, inter-individual heterogeneity is clearly crucial. Therefore, we introduce our novel approach, Standing Variations Modeling, which exploits inestimability of model parameters to capture heterogeneity in a deterministic model. We use data from mouse xenograft experiments to infer distributions on parameters critical to tumor growth and to the resultant immune response. Sampling model parameters from these distributions allows us to represent heterogeneity, both at the level of the tumor cells and the individual (mouse) being treated. Our model simulations explain the limited success of sipuleucel-T observed in practice, and predict an optimal combination regime that maximizes predicted efficacy. This approach will generalize to a range of emerging cancer immunotherapies.

A validated mathematical model of FGFR3-mediated tumor growth reveals pathways to harness the benefits of combination targeted therapy and immunotherapy in bladder cancer.

Bladder cancer is a common malignancy with over 80,000 estimated new cases and nearly 18,000 deaths per year in the United States alone. Therapeutic options for metastatic bladder cancer had not evolved much for nearly four decades, until recently, when five immune checkpoint inhibitors were approved by the U.S. Food and Drug Administration (FDA). Despite the activity of these drugs in some patients, the objective response rate for each is less than 25%. At the same time, fibroblast growth factor receptors (FGFRs) have been attractive drug targets for a variety of cancers, and in 2019 the FDA approved the first therapy targeted against FGFR3 for bladder cancer. Given the excitement around these new receptor tyrosine kinase and immune checkpoint targeted strategies, and the challenges they each may face on their own, emerging data suggest that combining these treatment options could lead to improved therapeutic outcomes. In this paper, we develop a mathematical model for FGFR3-mediated tumor growth and use it to investigate the impact of the combined administration of a small molecule inhibitor of FGFR3 and a monoclonal antibody against the PD-1/PD-L1 immune checkpoint. The model is carefully calibrated and validated with experimental data before survival benefits, and dosing schedules are explored. Predictions of the model suggest that FGFR3 mutation reduces the effectiveness of anti-PD-L1 therapy, that there are regions of parameter space where each monotherapy can outperform the other, and that pretreatment with anti-PD-L1 therapy always results in greater tumor reduction even when anti-FGFR3 therapy is the more effective monotherapy.

Receptor cross-talk in angiogenesis: mapping environmental cues to cell phenotype using a stochastic, Boolean signaling network model.

Cancer invasion and metastasis depend on tumor-induced angiogenesis, the means by which cancer cells attract and maintain a blood supply. During angiogenesis, cellular processes are tightly coordinated by signaling molecules and their receptors. Understanding how endothelial cells synthesize multiple biochemical signals can catalyze the development of novel therapeutic strategies to combat cancer. This study is the first to propose a signal transduction model highlighting the cross-talk between key receptors involved in angiogenesis, namely the VEGF, integrin, and cadherin receptors. From experimental data, we construct a network model of receptor cross-talk and analyze its dynamics. We identify relationships between receptor activation combinations and cellular function, and show that cross-talk is crucial to phenotype determination. The network converges to a unique set of output states that correspond to known cell phenotypes: migratory, proliferating, quiescent, apoptotic, and it predicts one phenotype that challenges the "go or grow" hypothesis. Finally, we use the model to study protein inhibition and to suggest molecular targets for anti-angiogenic therapies.

Topography of extracellular matrix mediates vascular morphogenesis and migration speeds in angiogenesis.

The extracellular matrix plays a critical role in orchestrating the events necessary for wound healing, muscle repair, morphogenesis, new blood vessel growth, and cancer invasion. In this study, we investigate the influence of extracellular matrix topography on the coordination of multi-cellular interactions in the context of angiogenesis. To do this, we validate our spatio-temporal mathematical model of angiogenesis against empirical data, and within this framework, we vary the density of the matrix fibers to simulate different tissue environments and to explore the possibility of manipulating the extracellular matrix to achieve pro- and anti-angiogenic effects. The model predicts specific ranges of matrix fiber densities that maximize sprout extension speed, induce branching, or interrupt normal angiogenesis, which are independently confirmed by experiment. We then explore matrix fiber alignment as a key factor contributing to peak sprout velocities and in mediating cell shape and orientation. We also quantify the effects of proteolytic matrix degradation by the tip cell on sprout velocity and demonstrate that degradation promotes sprout growth at high matrix densities, but has an inhibitory effect at lower densities. Our results are discussed in the context of ECM targeted pro- and anti-angiogenic therapies that can be tested empirically.