I will survey recent work in my group on multiscale models of cells. I will describe how we model the individual states of a cell, and what new features are introduced when cells interact, whether by chemical signaling, by mechanical forces, or by cell-cell adhesion. In some cases, the examples are taken from developmental biology, and in other cases from synthetic biology, wound healing, and metastasis. While my focus will be on computational models, some connections to nonlocal and continuum models will be made.

Mathematics underpins fundamental theories in physics such as quantum mechanics, general relativity, and quantum field theory. Nonetheless, its success in modern biology, namely cellular biology, molecular biology, chemical biology, genomics, and genetics, has been quite limited. Artificial intelligence (AI) has fundamentally changed the landscape of science, engineering, and technology in the past decade and holds a great future for discovering the rules of life. However, AI-based biological discovery encounters challenges arising from the intricate complexity, high dimensionality, nonlinearity, and multiscale biological systems. We tackle these challenges by a mathematical AI paradigm. We have introduced persistent cohomology, persistent spectral graphs, persistent path Laplacians, persistent sheaf Laplacians, and evolutionary de Rham-Hodge theory to significantly enhance AI's ability to tackle biological challenges. Using our mathematical AI approaches, my team has been the top winner in D3R Grand Challenges, a worldwide annual competition series in computer-aided drug design and discovery for years. By further integrating mathematical AI with millions of genomes isolated from patients, we uncovered the mechanisms of SARS-CoV-2 evolution and accurately forecast emerging dominant SARS-CoV-2 variants.

Stem cells are required for tissue maintenance and homeostasis during an organism's lifetime. When dividing, stem cells can either self-renew into stem cells, or their progeny can become progenitor cells that can then differentiate into more specialized cells. Feedback from the differentiated cell population onto regulation of division controls tissue growth and maintains tissue homeostasis. Here I consider how to differentiate between multiple cell lineage models with potential nonlinear feedback terms. I consider the influence of several commonly made assumptions in stem cell models including: (1) whether or not more differentiated progeny can divide, (2) whether or not stem cell death is included, (3) the impact of symmetric and asymmetric stem cell divisions. Including differentiated cell division in the simplified model can lead to emergence of a spurious steady state that may not be biologically realistic, and the parameter region for existence of nontrivial equilibrium, which corresponds to tissue at homeostasis, shrinks rapidly as more feedbacks are added to the model. I consider two potential ways to modify a stem cell lineage model to get rid of this unrealistic steady state: getting rid of differentiated cell division altogether or modifying their division to depend on stem cells.
Next, I consider a model for neural stem cells (NSCs) that can be in an actively dividing state or in a quiescent state. The balance between stem cell quiescence and cycling activity determines the rate of neurogenesis. With age, more NSCs enter the quiescent state, while the total number of NSCs decreases. I consider which model variant could best explain the observed decline of neural stem cells seen in experimental data from mice.

The creation of digital twins (DTs) of cancer patients can assist us in predicting the evolution of an individual's cancer through modeling each tumor’s characteristics and response to treatment. We therefore take advantage of new advances in computational approaches and combine mechanistic, machine learning, and stochastic modeling approaches to create a DT platform, which utilizes biological, biomedical, and EHR data sets. For each patient, the DT receives their information as input and predicts the evolution of their cancer. We propose to develop a mechanistic model based on the quantitative systems pharmacology (QSP) modeling, which is one of the main computational approaches used to discover, test, and predict dose-exposure response. One of the main challenges of the QSP modeling is parameter estimation. Traditionally, these models assume all patients have similar diseases, and the values of parameters of the QSP model are identical for all patients. Therefore, parameters are commonly calibrated using the data often assembled from disparate sources. To develop a personalized DT, we use patient-specific data for parameter estimations, sensitivity analysis, and uncertainty quantification. For each patient, we estimate the values of parameters of their QSP model using their data. We perform a multi-dimensional sensitivity analysis and uncertainty quantification on the mechanistic model to find a set of critical interactions and predict the intervals of confidence. Since this QSP model includes the data-driven mechanistic model of cells and molecules' interaction networks, one of the ultimate results of this DT is the prediction of evolution of cancer in response to a given targeted therapy.

Drug discovery is a molecular search task with a complex objective: modify the function of a complex biological system to interrupt disease processes. Conventionally, it is a costly, high failure-rate process – with molecular candidates clearing preclinical safety and efficacy hurdles only to fail upon delivery to humans. This happens because early screens in the discovery pipeline fall short of capturing the ultimate therapeutic value of new molecular candidates. For a molecule to become a successful drug, it should: (1) bind to a desired target protein; (2) be deliverable from a desired site of administration (oral, intravenous, etc.) to the physiological site of activity, with sufficient concentration for a sufficient duration of time; and (3) promote the desired pharmacological effect without causing unwanted toxicity. The chemical space that meets one of these objectives likely requires compromises in another, as binding, delivery, and activity depend on coupled and dynamic biophysical and biochemical interactions. Thus, to improve the success rate of drug discovery, we should ideally look at design through the lens of human physiology. Quantitative systems pharmacology models could offer the molecule-to-therapeutic-outcome mapping required to inform AI-driven drug design. However, these models present multiple challenges – from the complexity of the biological pathways driving disease processes to the knowledge gaps limiting model construction and parameterization, to the challenges presented by data limitations and the relative computational expense of mechanistic systems models. In this talk, I will present our work on enabling physiology-informed, AI-driven design of new molecular entities.

Neurodegenerative diseases such as Alzheimer’s or Parkinson’s are devastating conditions with poorly understood mechanisms and no known cure. Yet a striking feature of these conditions is the characteristic pattern of invasion throughout the brain, leading to well-codified disease stages visible to neuropathology and associated with various cognitive deficits and pathologies. This evolution is associated with the aggregation of key toxic proteins. In this talk, I will show how we use multiscale modelling to gain insight into this process and, doing so, gain some understanding on how the brain works. In particular, by looking at protein dynamics on the neuronal network, we we can unravel some of the universal features associated with dementia that are driven by both the network topology and protein kinetics. By further coupling this approach with functional models of the brain, we will show that we can explain important aspects of cognitive loss during disease development.

In this talk I will present a general theory for active fluids which convert chemical energy into various type of mechanical energy. This is the extension of the classical energetic variational approaches for mechanical systems. The methods will cover a range of both chemical reaction kenetics and mechanical processes.

Alzheimer's disease (AD) affects over 5 million people in the US, and personalized treatment plans have emerged as a promising approach to managing the disease. However, analyzing the growing electronic AD brain data requires a new approach. This talk introduces a mathematical modeling approach that describes the progression of AD clinical biomarkers and incorporates patient data for personalized prediction and treatment optimization. The mathematical model is validated on a multi-institutional dataset of AD biomarkers, such as the Alzheimer's Disease Neuroimaging Initiative (ADNI) dataset, and personalized therapeutic simulation studies are performed using optimal controls. The results of this study can provide valuable insights for personalized AD treatment planning.

Discovering in the 1830's the predominance of Fibonacci numbers in the number of spirals in plants might have been the first true bio-mathematical work in history. Plants that do not fit that pattern have often be classified as exceptions, belonging to other Fibonacci-like sequences, or just plain ignored. It turns out that many of these patterns can be classified into a general category that our group has called Quasi-symmetric, that exhibit numbers of spirals of the form varying between (n, n), (n, n+1) and sometimes (n, n+2) (e.g. strawberries, corn and many more). In this talk I will show how one can reproduce the dichotomy between Fibonacci and Quasi-symmetric by tuning of a growth parameter in a simple computational dynamical model. I will also show how one can systematically collect and analyze plant data with that dichotomy in mind.

The brain tissue is protected by the blood-brain barrier: a wall of tightly-packed cells that keep unwanted molecules from crossing from the blood vessels into the tissue. This presents challenges to delivering therapeutic drugs to locations in the brain to treat certain diseases. One approach to meeting this challenge is to encapsulate the drugs in sono-sensitive nano-carriers. These vesicles can then be made to release their cargo locally using focused ultra-sound beams at intensities that are not damaging to the surrounding tissue. Mathematical problems come up when trying to answer questions such as: How can we describe the kinetics of drug transport and distribution through the tissue? What is the best positioning for an array of ultrasound transducers in order to produce the required signal at the right spot in the brain? What ultrasound parameters and dosages produce the the desired drug profile at the target region? 

In this talk I will discuss the specific mathematical challenges, as well as some approaches to their solution. 


In the first part of the talk, I will present an array of microfluidic and 3D culture platforms that can be used to study the effects of microenvironmental factors on tumor metastasis with a focus on macrophages. I will show how tumor-derived paracrine signals and extracellular matrix remodeling play a critical role in recruiting macrophages that in turn establish a pro-invasive and drug-resistant microenvironment. In the second part of the talk, I will present an integrated experimental-computational approach to study the effects of fibroblasts on HER2-targeted therapy drug resistance in breast cancer. A two-cell circuit model was developed using dynamic measurements of tumor cell growth and death rates. To identify molecular mechanisms associated with fibroblast-mediated drug resistance we are profiling signaling pathway activity using multiplexed in situ immunofluorescence staining and reverse phase proteomics. The results from the proposed studies will lead to the identification of microenvironmental niches and associated tumor cell signaling pathways that could serve as biomarkers for patient stratification and provide important information for the design of rational combination therapies that will re-sensitize tumors to treatment.

Malignant gliomas are highly invasive brain tumors. Recent attention has focused on their capacity for network-driven invasion, whereby mitotic events can be followed by the migration of nuclei along long thin cellular protrusions, termed tumour microtubes (TM). Here I develop a mathematical model that describes this microtube-driven invasion of gliomas. I show that scaling limits lead to well known glioma models as special cases such as go-or-grow models, the PI model of Swanson, and the anisotropic model of Swan. Numerical simulations are used to compare between the models. (Joint work with N. Loy, K.J. Painter, R. Thiessen).

We will discuss approaches we have taken to mathematically modeling immune dynamics in the context of cancer, type I diabetes, and viral infection (SARS-CoV2).

The difficulty of predicting the emergent development, homeostasis and disfunction of tissues from cells’ molecular signatures limits our ability to integrate molecular and genetic information to make meaningful predictions at the organ or organism level. Virtual Tissues are an approach to constructing quantitative, predictive Agent-Based mechanistic models starting from cell behaviors and combining subcellular molecular kinetics models, the physical and mechanical behaviors of cells and the longer-range effects of the extracellular environment. For the past 15 years, we have been developing the open-source Virtual-Tissue model specification and execution framework CompuCell3D (CC3D) (www.compucell3d.org) which aims to make Virtual-Tissue modeling more accessible to biologists and bioengineers. I will talk about some recent applications of CC3D(to model in host infection and immune response and corneal damage and recovery and some extensions of CC3D to function better with the increasingly popular Jupyter Notebook environment. CC3D uses the lattice-based Cellular-Potts Model (CPM) formalism to represent cells, which allows for description of cell shape but also has certain limitations. Recently, we have developed a scriptable Python-native (Jupyter-Notebook) off-lattice Center and Vertex model Virtual-Tissue modeling environment, Tissue Forge (https://tissue-forge-documentation.readthedocs.io/en/latest/) which allows for real-time interactive simulations with hundreds of thousands of cells. I will introduce Tissue Forge with the example of modeling zebrafish epiboly. Finally, solving the diffusion equation is one of the most computationally costly components of many Virtual Tissue models. I will discuss our preliminary efforts to develop U-net based Machine-Learning (ML) surrogates for diffusion solvers and invite discussion of how and where to integrate ML and mechanistic modeling.

A central feature of circadian systems is their response to an external, pacemaking 24 hour light-dark (LD) drive which typically leads to entrainment of circadian oscillator. There are, however, several naturally arising situations in which a circadian system is incapable of entrainment, either due to abnormal intrinsic properties of the oscillators or due to changes in the light-dark input that the oscillators receive. In this talk, we will use entrainment maps to describe circumstances that fall outside the normal fixed phase relationship between LD forcing and oscillator such as during jet lag, shift work and non-24 hour sleep-wake disorder. The mathematical and computational methods used to study these problems revolve around finding stable limit cycle solutions of the governing equations and it is the reduction of this study to a one-dimensional framework that will be the focus of the talk.

Small GTPases, such as Rac and Rho, are well known central regulators of cell morphology and motility, whose dynamics also play a role in coordinating collective cell migration. Experiments have shown GTPase dynamics to be affected by both chemical and mechanical cues, but also to be spatially and temporally heterogeneous. This heterogeneity is found both within a single cell, and between cells in a tissue. For example, sometimes the leader and follower cells display an inverted GTPase configuration. While progress on understanding GTPase dynamics in single cells has been made, a major remaining challenge is to understand the role of GTPase heterogeneity in collective cell migration. Motivated by recent one-dimensional experiments (e.g. micro-channels) we introduce a one-dimensional modelling framework allowing us to integrate cell bio-mechanics, changes in cell size, and detailed intra-cellular signalling circuits (reaction-diffusion equations). Using this framework, we build cell migration models of both loose (mesenchymal) and cohering (epithelial) tissues. We use numerical simulations, and analysis tools, such as local perturbation analysis, to provide insights into the regulatory mechanisms coordinating collective cell migration. We show how feedback from mechanical tension to GTPase activation lead to a variety of dynamics, resembling both normal and pathological behavior.

Actin cytoskeleton networks generate local topological signatures due to the natural variations in the number, size, and shape of holes of the networks. Persistent homology is a method that explores these topological properties of data and summarizes them as persistence diagrams. In this work, we analyze and classify these filament networks by transforming them into persistence diagrams whose variability is quantified via a Bayesian framework on the space of persistence diagrams. The proposed generalized Bayesian framework adopts an independent and identically distributed cluster point process characterization of persistence diagrams and relies on a substitution likelihood argument. This framework provides the flexibility to estimate the posterior cardinality distribution of points in a persistence diagram and the posterior spatial distribution simultaneously. We present a closed form of the posteriors under the assumption of Gaussian mixtures and binomials for prior intensity and cardinality respectively. Using this posterior calculation, we implement a Bayes factor algorithm to classify the actin filament networks and benchmark it against several state-of-the-art classification methods.

Astrocytes are glial cells making up 50% of brain volume, and playing multiple important roles. We are developing tools to include “effective” astrocytes in neuronal network models in an easy-to-implement, and relatively computationally-efficient way. I will show results from modeling different aspects of astrocytes at several vastly different spatial and temporal scales. We will discuss the effects of desensitization of astrocyte G-protein-coupled receptors, neurotransmitter diffusion with recharging traps (DiRT), and modulation of neuronal excitability due to changes in extracellular ion concentrations. From these detailed models we extract the essential ways in which astrocytes influence nearby neurons, and include these effects in developing an updated understanding of neural network dynamics and synchronization.

Tumors are typically comprised of heterogeneous cell populations exhibiting diverse phenotypes. This heterogeneity, which is correlated with tumor aggressiveness and treatment-failure, confounds current drug screening efforts to identify effective candidate therapies for individual tumors. In the first part of the talk I will present a modeling-driven statistical framework that enables the deconvolution of tumor samples into individual subcomponents exhibiting differential drug-response, using standard bulk drug-screen measurements. In the second part of the talk I will present some efforts towards obtaining insights about tumor evolution from standard genomic data. In particular, we analyze the site frequency spectrum (SFS), a population summary statistic of genomic data, for exponentially growing tumor populations, and we demonstrate how these results can in principle be used to gain insights into tumor evolutionary parameters.

All organisms are composed of multiple chemical elements such as nitrogen (N), phosphorus (P), and carbon (C). P is essential to build nucleic acids (DNA and RNA) and N is needed for protein production. To keep track of the mismatch between P requirement in the consumer and the P content in the producer, stoichiometric models have been constructed to explicitly incorporate food quality and quantity. Most stoichiometric models have suggested that the consumer dynamics heavily depend on P content in the producer when the producer has low nutrient content (low P:C ratio). Motivated by recent lab experiments, researchers explored the effect of excess producer nutrient content (extremely high P:C ratio) on the grazer dynamics. This phenomenon is called the stoichiometric knife edge. However, the global analysis of these resource quality based models is challenging because the phase plane/space is separated into many regions in which the governing nonlinear equations are different. The aim of this talk is to present an overview of the rich and novel dynamics embodied in these stoichiometric population models and its many biological implications and present an alternative framework to build mathematically more tractable and biologically more plausible models.

In this presentation, I will discuss the construction and results of numerical simulations quantifying flows around several species of soft corals. In the first project, the flows near the tentacles of xeniid soft corals are quantified for the first time. Their active pulsations are thought to enhance their symbionts' photosynthetic rates by up to an order of magnitude. These polyps are approximately 1 cm in diameter and pulse at frequencies between approximately 0.5 and 1 Hz. As a result, the frequency-based Reynolds number calculated using the tentacle length and pulse frequency is on the order of 10 and rapidly decays as with distance from the polyp. This introduces the question of how these corals minimize the reversibility of the flow and bring in new volumes of fluid during each pulse. We estimate the Péclet number of the bulk flow generated by the coral as being on the order of 100–1000 whereas the flow between the bristles of the tentacles is on the order of 10. This illustrates the importance of advective transport in removing oxygen waste. In the second project, the flows through the elaborate branching structures of gorgonian colonies are considered. As water moves through the elaborate branches, it is slowed, and recirculation zones can form downstream of the colony. At the smaller scale, individual polyps that emerge from the branches expand their tentacles, further slowing the flow. At the smallest scale, the tentacles are covered in tiny pinnules where exchange occurs. We quantified the gap to diameter ratios for various gorgonians at the scale of the branches, the polyp tentacles and the pinnules. We then used computational fluid dynamics to determine the flow patterns at all three levels of branching. We quantified the leakiness between the branches, tentacles and pinnules over the biologically relevant range of Reynolds numbers and gap-to-diameter ratios and found that the branches and tentacles can act as either leaky rakes or solid plates depending upon these dimensionless parameters. The pinnules, in contrast, mostly impede the flow. Using an agent-based modeling framework, we quantified plankton capture as a function of the gap-to diameter ratio of the branches and the Reynolds number. We found that the capture rate depends critically on both morphology and Reynolds number.