Understanding the intricate interactions of cancer cells with the tumor microenvironment (TME) is a pre-requisite for the optimization of immunotherapy. Mechanistic models such as quantitative systems pharmacology (QSP) provide insights into the TME dynamics and predict the efficacy of immunotherapy in virtual patient populations/digital twins but require vast amounts of multimodal data for parameterization. Large-scale datasets characterizing the TME are available due to recent advances in bioinformatics for multi-omics data. Here, we discuss the perspectives of leveraging omics-derived bioinformatics estimates to inform QSP models and circumvent the challenges of model calibration and validation in immuno-oncology.

Advancements in imaging technologies have revolutionized our ability to deeply profile pathological tissue architectures, generating large volumes of imaging data with unparalleled spatial resolution. This type of data collection, namely, spatial proteomics, offers invaluable insights into various human diseases. Simultaneously, computational algorithms have evolved to manage the increasing dimensionality of spatial proteomics inherent in this progress. Numerous imaging-based computational frameworks, such as computational pathology, have been proposed for research and clinical applications. However, the development of these fields demands diverse domain expertise, creating barriers to their integration and further application. This review seeks to bridge this divide by presenting a comprehensive guideline. We consolidate prevailing computational methods and outline a roadmap from image processing to data-driven, statistics-informed biomarker discovery. Additionally, we explore future perspectives as the field moves toward interfacing with other quantitative domains, holding significant promise for precision care in immuno-oncology.

Quantitative systems pharmacology (QSP) models and spatial agent-based models (ABM) are powerful and efficient approaches for the analysis of biological systems and for clinical applications. Although QSP models are becoming essential in discovering predictive biomarkers and developing combination therapies through in silico virtual trials, they are inadequate to capture the spatial heterogeneity and randomness that characterize complex biological systems, and specifically the tumor microenvironment. Here, we extend our recently developed spatial QSP (spQSP) model to analyze tumor growth dynamics and its response to immunotherapy at different spatio-temporal scales. In the model, the tumor spatial dynamics is governed by the ABM, coupled to the QSP model, which includes the following compartments: central (blood system), tumor, tumor-draining lymph node, and peripheral (the rest of the organs and tissues). A dynamic recruitment of T cells and myeloid-derived suppressor cells (MDSC) from the QSP central compartment has been implemented as a function of the spatial distribution of cancer cells. The proposed QSP-ABM coupling methodology enables the spQSP model to perform as a coarse-grained model at the whole-tumor scale and as an agent-based model at the regions of interest (ROIs) scale. Thus, we exploit the spQSP model potential to characterize tumor growth, identify T cell hotspots, and perform qualitative and quantitative descriptions of cell density profiles at the invasive front of the tumor. Additionally, we analyze the effects of immunotherapy at both whole-tumor and ROI scales under different tumor growth and immune response conditions. A digital pathology computational analysis of triple-negative breast cancer specimens is used as a guide for modeling the immuno-architecture of the invasive front.

Annotating cancerous regions in whole-slide images (WSIs) of pathology samples plays a critical role in clinical diagnosis, biomedical research, and machine learning algorithms development. However, generating exhaustive and accurate annotations is labor-intensive, challenging, and costly. Drawing only coarse and approximate annotations is a much easier task, less costly, and it alleviates pathologists' workload. In this paper, we study the problem of refining these approximate annotations in digital pathology to obtain more accurate ones. Some previous works have explored obtaining machine learning models from these inaccurate annotations, but few of them tackle the refinement problem where the mislabeled regions should be explicitly identified and corrected, and all of them require a - often very large - number of training samples. We present a method, named Label Cleaning Multiple Instance Learning (LC-MIL), to refine coarse annotations on a single WSI without the need for external training data. Patches cropped from a WSI with inaccurate labels are processed jointly within a multiple instance learning framework, mitigating their impact on the predictive model and refining the segmentation. Our experiments on a heterogeneous WSI set with breast cancer lymph node metastasis, liver cancer, and colorectal cancer samples show that LC-MIL significantly refines the coarse annotations, outperforming state-of-the-art alternatives, even while learning from a single slide. Moreover, we demonstrate how real annotations drawn by pathologists can be efficiently refined and improved by the proposed approach. All these results demonstrate that LC-MIL is a promising, lightweight tool to provide fine-grained annotations from coarsely annotated pathology sets.

Macrophages respond to signals in the microenvironment by changing their functional phenotypes, a process known as polarization. Depending on the context, they acquire different patterns of transcriptional activation, cytokine expression and cellular metabolism which collectively constitute a continuous spectrum of phenotypes, of which the two extremes are denoted as classical (M1) and alternative (M2) activation. To quantitatively decode the underlying principles governing macrophage phenotypic polarization and thereby harness its therapeutic potential in human diseases, a systems-level approach is needed given the multitude of signaling pathways and intracellular regulation involved. Here we develop the first mechanism-based, multi-pathway computational model that describes the integrated signal transduction and macrophage programming under M1 (IFN-γ), M2 (IL-4) and cell stress (hypoxia) stimulation. Our model was calibrated extensively against experimental data, and we mechanistically elucidated several signature feedbacks behind the M1-M2 antagonism and investigated the dynamical shaping of macrophage phenotypes within the M1-M2 spectrum. Model sensitivity analysis also revealed key molecular nodes and interactions as targets with potential therapeutic values for the pathophysiology of peripheral arterial disease and cancer. Through simulations that dynamically capture the signal integration and phenotypic marker expression in the differential macrophage polarization responses, our model provides an important computational basis toward a more quantitative and network-centric understanding of the complex physiology and versatile functions of macrophages in human diseases.

Vascular endothelial growth factor (VEGF) is a powerful regulator of neovascularization. VEGF binding to its cognate receptor, VEGFR2, activates a number of signaling pathways including ERK1/2. Activation of ERK1/2 is experimentally shown to involve sphingosine kinase 1 (SphK1) activation and its calcium-dependent translocation downstream of ERK1/2. Here we construct a rule-based computational model of signaling downstream of VEGFR2, by including SphK1 and calcium positive feedback mechanisms, and investigate their consequences on ERK1/2 activation. The model predicts the existence of VEGF threshold in ERK1/2 activation that can be continuously tuned by cellular concentrations of SphK1 and sphingosine 1 phosphate (S1P). The computer model also predicts powerful effects of perturbations in plasma and ER calcium pump rates and the current through the CRAC channels on ERK1/2 activation dynamics, highlighting the critical role of intracellular calcium in shaping the pERK1/2 signal. The model is then utilized to simulate anti-angiogenic therapeutic interventions targeting VEGFR2-ERK1/2 axis. Simulations indicate that monotherapies that exclusively target VEGFR2 phosphorylation, VEGF, or VEGFR2 are ineffective in shutting down signaling to ERK1/2. By simulating therapeutic strategies that target multiple nodes of the pathway such as Raf and SphK1, we conclude that combination therapy should be much more effective in blocking VEGF signaling to EKR1/2. The model has important implications for interventions that target signaling pathways in angiogenesis relevant to cancer, vascular diseases, and wound healing.

Hypoxia is an important physiological stress signal that drives angiogenesis, the formation of new blood vessels. Besides an increase in the production of pro-angiogenic signals such as vascular endothelial growth factor (VEGF), hypoxia also stimulates the production of anti-angiogenic signals. Thrombospondin-1 (TSP-1) is one of the anti-angiogenic factors whose synthesis is driven by hypoxia. Cellular synthesis of TSP-1 is tightly regulated by different intermediate biomolecules including proteins that interact with hypoxia-inducible factors (HIFs), transcription factors that are activated by receptor and intracellular signaling, and microRNAs which are small non-coding RNA molecules that function in post-transcriptional modification of gene expression. Here we present a computational model that describes the mechanistic interactions between intracellular biomolecules and cooperation between signaling pathways that together make up the complex network of TSP-1 regulation both at the transcriptional and post-transcriptional level. Assisted by the model, we conduct in silico experiments to compare the efficacy of different therapeutic strategies designed to modulate TSP-1 synthesis in conditions that simulate tumor and peripheral arterial disease microenvironment. We conclude that TSP-1 production in endothelial cells depends on not only the availability of certain growth factors but also the fine-tuned signaling cascades that are initiated by hypoxia.

HRMs (hypoxia-responsive miRNAs) are a specific group of microRNAs that are regulated by hypoxia. Recent studies revealed that several HRMs including let-7 family miRNAs were highly induced in response to HIF (hypoxia-inducible factor) stabilization in hypoxia, and they potently participated in angiogenesis by targeting AGO1 (argonaute 1) and upregulating VEGF (vascular endothelial growth factor). Here we constructed a novel computational model of microRNA control of HIF-VEGF pathway in endothelial cells to quantitatively investigate the role of HRMs in modulating the cellular adaptation to hypoxia. The model parameters were optimized and the simulations based on these parameters were validated against several published in vitro experimental data. To advance the mechanistic understanding of oxygen sensing in hypoxia, we demonstrated that the rate of HIF-1α nuclear import substantially influences its stabilization and the formation of HIF-1 transcription factor complex. We described the biological feedback loops involving let-7 and AGO1 in which the impact of external perturbations were minimized; as a pair of master regulators when low oxygen tension was sensed, they coordinated the critical process of VEGF desuppression in a controlled manner. Prompted by the model-motivated discoveries, we proposed and assessed novel pathway-specific therapeutics that modulate angiogenesis by adjusting VEGF synthesis in tumor and ischemic cardiovascular disease. Through simulations that capture the complex interactions between miRNAs and miRNA-processing molecules, this model explores an innovative perspective about the distinctive yet integrated roles of different miRNAs in angiogenesis, and it will help future research to elucidate the dysregulated miRNA profiles found in cancer and various cardiovascular diseases.

Vascular endothelial growth factor (VEGF) is a potent cytokine that binds to specific receptors on the endothelial cells lining blood vessels. The signaling cascade triggered eventually leads to the formation of new capillaries, a process called angiogenesis. Distributions of VEGF receptors and VEGF ligands are therefore crucial determinants of angiogenic events and, to our knowledge, no quantification of abluminal vs. luminal receptors has been performed. We formulate a molecular-based compartment model to investigate the VEGF distribution in blood and tissue in humans and show that such quantification would lead to new insights on angiogenesis and VEGF-dependent diseases. Our multiscale model includes two major isoforms of VEGF (VEGF(121) and VEGF(165)), as well as their receptors (VEGFR1 and VEGFR2) and the non-signaling co-receptor neuropilin-1 (NRP1). VEGF can be transported between tissue and blood via transendothelial permeability and the lymphatics. VEGF receptors are located on both the luminal and abluminal sides of the endothelial cells. In this study, we analyze the effects of the VEGF receptor localization on the endothelial cells as well as of the lymphatic transport. We show that the VEGF distribution is affected by the luminal receptor density. We predict that the receptor signaling occurs mostly on the abluminal endothelial surface, assuming that VEGF is secreted by parenchymal cells. However, for a low abluminal but high luminal receptor density, VEGF binds predominantly to VEGFR1 on the abluminal surface and VEGFR2 on the luminal surface. Such findings would be pertinent to pathological conditions and therapies related to VEGF receptor imbalance and overexpression on the endothelial cells and will hopefully encourage experimental receptor quantification for both luminal and abluminal surfaces on endothelial cells.

Members of the vascular endothelial growth factor (VEGF) family of proteins are critical regulators of angiogenesis. VEGF concentration gradients are important for activation and chemotactic guidance of capillary sprouting, but measurement of these gradients in vivo is not currently possible. We have constructed a biophysically and molecularly detailed computational model to study microenvironmental transport of two isoforms of VEGF in rat extensor digitorum longus skeletal muscle under in vivo conditions. Using parameters based on experimental measurements, the model includes: VEGF secretion from muscle fibers; binding to the extracellular matrix; binding to and activation of endothelial cell surface VEGF receptors; and internalization. For 2-D cross sections of tissue, we analyzed predicted VEGF distributions, gradients, and receptor binding. Significant VEGF gradients (up to 12% change in VEGF concentration over 10 mum) were predicted in resting skeletal muscle with uniform VEGF secretion, due to non-uniform capillary distribution. These relative VEGF gradients were not sensitive to extracellular matrix composition, or to the overall VEGF expression level, but were dependent on VEGF receptor density and affinity, and internalization rate parameters. VEGF upregulation in a subset of fibers increased VEGF gradients, simulating transplantation of pro-angiogenic myoblasts, a possible therapy for ischemic diseases. The number and relative position of overexpressing fibers determined the VEGF gradients and distribution of VEGF receptor activation. With total VEGF expression level in the tissue unchanged, concentrating overexpression into a small number of adjacent fibers can increase the number of capillaries activated. The VEGF concentration gradients predicted for resting muscle (average 3% VEGF/10 mum) is sufficient for cellular sensing; the tip cell of a vessel sprout is approximately 50 mum long. The VEGF gradients also result in heterogeneity in the activation of blood vessel VEGF receptors. This first model of VEGF tissue transport and heterogeneity provides a platform for the design and evaluation of therapeutic approaches.

Angiogenesis (neovascularization) plays a crucial role in a variety of physiological and pathological conditions including cancer, cardiovascular disease, and wound healing. Vascular endothelial growth factor (VEGF) is a critical regulator of angiogenesis. Multiple VEGF receptors are expressed on endothelial cells, including signaling receptor tyrosine kinases (VEGFR1 and VEGFR2) and the nonsignaling co-receptor Neuropilin-1. Neuropilin-1 binds only the isoform of VEGF responsible for pathological angiogenesis (VEGF165), and is thus a potential target for inhibiting VEGF signaling. Using the first molecularly detailed computational model of VEGF and its receptors, we have shown previously that the VEGFR-Neuropilin interactions explain the observed differential effects of VEGF isoforms on VEGF signaling in vitro, and demonstrated potent VEGF inhibition by an antibody to Neuropilin-1 that does not block ligand binding but blocks subsequent receptor coupling. In the present study, we extend that computational model to simulation of in vivo VEGF transport and binding, and predict the in vivo efficacy of several Neuropilin-targeted therapies in inhibiting VEGF signaling: (a) blocking Neuropilin-1 expression; (b) blocking VEGF binding to Neuropilin-1; (c) blocking Neuropilin-VEGFR coupling. The model predicts that blockade of Neuropilin-VEGFR coupling is significantly more effective than other approaches in decreasing VEGF-VEGFR2 signaling. In addition, tumor types with different receptor expression levels respond differently to each of these treatments. In designing human therapeutics, the mechanism of attacking the target plays a significant role in the outcome: of the strategies tested here, drugs with similar properties to the Neuropilin-1 antibody are predicted to be most effective. The tumor type and the microenvironment of the target tissue are also significant in determining therapeutic efficacy of each of the treatments studied.