## 2016-2017 Seminars

### Postdoctoral Seminars

The world of genetics is highly amenable to statistical analysis. Genetical processes are random; they tend to have clear-cut outcomes that are few in quantity; and they almost always occur in vast numbers. In this talk I'll discuss some of the ways in which statistics can help us understand the incomprehensibly tiny, but nevertheless vitally important, world of genetics.

Understanding the cell dynamics is important to determine the origin of many diseases including cancer. It can also help us to improve treatments by suggesting ways to control the growth rate of cells and minimizing the spread of mutants. In this talk, we introduce and compare available computational models for cell dynamics in colon and intestinal crypts. Then, we introduce a new stochastic model for colon and intestinal crypt, which incorporates the new experimental discoveries. Using this model, we study the spread of mutants in colon and intestinal crypt. The model suggests ways to improve the outcome of stem cell therapies.

Normal mode analysis, also known as harmonic analysis, is applicable to many fields. It identifies the natural, resonant movements of a physical object such as a building, guitar string, or molecule. In this talk I will demonstrate how to apply normal mode analysis to Gag proteins, the building block of HIV. Specifically, this technique will be used to analyze how the movement of the proteins change when bound to specific RNA during HIV viral assembly.

The spread of disease in a household can be represented as a SID model (susceptible-infective-disease). We focus on disease propagation at the household level with age-structure and environment contamination. Using cross-sectional data, we estimate parameters for a likelihood representative of the household structure. Using this likelihood, we generate data that can be used to fit the SID model in steady state via MCMC methods. This provides a template for fitting dynamic models using static information.

In this talk, we propose a new approach for clustering objects in general datasets without constraints on the dimensions. It combines ideas from principal component analysis (PCA), outlier filtering, and von Mises-Fisher mixture models. The proposed algorithm can also help select the optimal number of clusters for objects of interest. Furthermore, we review several methods on choosing the number of groups in clustering, including the Gap statistics and average silhouette method. An extensive Monte Carlo simulation study is also performed to compare the proposed approach with other popular ones on estimating the number of clusters. Finally, we apply this novel approach to the omics datasets involving the genes in interleukin-4 (IL4) signaling pathway for various types of cancer, and interesting results are obtained which are consistent with previous findings in literature.

Mass cytometers can record tens of features for millions of cells in a sample, and in particular, for leukemic cells. Many methods consider how to cluster or identify populations of phenotypically similar cells within cytometry data, but there has yet to be a connection between cell activity and other features and these groups or clusters. We use differential geometric ideas to consider how cell cycle and signaling features vary as a function of the cell populations. This consideration leads to a better understanding of the nonlinear relationships that exist in the cytometry data.

Lung cancer, primarily non-small-cell lung cancer (NSCLC), is the leading cause of cancer deaths in the United States and worldwide. While early detection significantly improves five-year survival, there are no reliable diagnostic tools for early detection. Several exosomal microRNAs (miRs) are overexpressed in NSCLC, and have been suggested as potential biomarkers for early detection. In this talk, I want to talk about a mathematical model for early stage of NSCLC with emphasis on the role of the three highest overexpressed miRs, namely miR-21, miR-205 and miR-155. Simulations of the model provide quantitative relationships between the tumor volume and the total mass of each of the above miRs in the tumor. Because of the positive correlation between these miRs in the tumor tissue and in the blood, the results may be viewed as a step toward establishing miRs 21, 205 and 155 as reliable serum biomarkers for early detection of NSCLC.

Phylogenetic statistical models can be used to determine evolutionary relationships between species. These models have the structure of algebraic varieties and understanding their geometry is important for proving many of their properties. For various biological reasons, phylogenetic *mixture models* are often a better model for species evolution. The variety of probability distributions associated to an r-tree mixture model is in general a *join variety*. In this talk, we examine some of the algebraic geometry underlying phylogenetic mixture models. In particular, we will show some computational tools for determining their generating sets.

The transmission of avian influenza between humans is extremely rare, and it mostly affects individuals who are in contact with infected family member. Although this scenario is uncommon, there have been multiple outbreaks that occur in small infection clusters with relatively low transmissibility, and thus are too weak to cause an epidemic. Still, subcritical transmission from stuttering chain data is vital for determining whether avian influenza is close to the threshold of R0 >1. In this presentation, we will explore two methods of estimating R0 using transmission chains and parameter estimation through data fitting.

The ability to produce force is fundamental to life at every scale. A great deal is understood about the individual molecules involved in force production, but far less is known about how they coordinate to produce whole-cell-level effects. In this talk, I present a simple experimental setup for producing cells with tunable geometric properties, and explore a mathematical model for how variations in cell geometry lead to changes in cell mechanical behavior, including the magnitude of force produced and the stability---or lack of stability---of cell symmetry.

Banded vegetation patterns have been observed in water-limited ecosystems across the globe, often appearing along sloped terrain with the stripes aligned transverse to the grade. Initial field observations suggest a link between the arcing of the bands and the curvature in elevation transverse to the slope. We present a simple model that qualitatively captures the influence of topography on the vegetation pattern. In particular, the model predicts that the curvature of the stripes can change direction when growing on top of a ridge versus within a channel. The model also predicts transitions between qualitatively different striped patterns as a function of the channel depth, precipitation parameter and plant mortality parameter.

Neural networks generate a variety of rhythmic activity patterns, often involving different timescales. One example arises in the respiratory network in the preBötzinger complex of the mammalian brainstem, which can generate the eupneic rhythm associated with normal respiration as well as recurrent low-frequency, large-amplitude bursts associated with sighing. Two competing hypotheses have been proposed to explain sigh generation: the recruitment of a neuronal population distinct from the eupneic rhythm-generating subpopulation or the reconfiguration of activity within a single population. Here, we consider two recent computational models, one of which represents each of the hypotheses. We use methods of dynamical systems theory, such as fast-slow decomposition, averaging, and bifurcation analysis, to understand the multiple timescale mechanisms underlying sigh generation in each model. In the course of our analysis, we discover that a third timescale is required to generate sighs in both models. Furthermore, we identify the similarities of the underlying mechanisms in the two models and the aspects in which they differ.

The nuclear pore complex controls all transport between the nucleus and the cytoplasm. Most steps of this transport process are driven by thermal fluctuations, not by an active input of energy, and yet the nuclear pore complex achieves high throughput of cargo while efficiently filtering out macromolecules that should not be transported. In this talk, I will outline what is known and what is still a mystery about how the nuclear pore functions, and I will discuss my preliminary efforts to model the mechanics of transport through the pore.

### Visitor Seminars

I will describe some ideas related to Persistent Homology by following an example related to the structure of the hippocampal neural code.

Base-calling is the signal processing step to identify the DNA sequences based on the raw signals measured by high-throughput DNA Sequencing machine. The accuracy of base-calling is crucial for high-throughput DNA sequencing and downstream analysis. Accordingly, we made an endeavor to reduce DNA sequencing errors of Illumina systems by correcting three kinds of crosstalk in the cluster intensity data. We discovered that signal crosstalk between adjacent clusters accounts for a large portion of sequencing errors in Illumina systems, even after correcting color crosstalk caused by the overlap of dye emission spectra and phasing/pre-phasing caused by out-of-step nucleotide synthesis. Interestingly and importantly, the spatial crosstalk between adjacent clusters is cluster-specific and often asymmetric, which cannot be corrected by existing deconvolution methods. Therefore, we introduce a novel mathematical method able to estimate and remove spatial crosstalk, thereby reducing base-calling errors by 44-69%. Furthermore, the resolution gained from this study provides new room for higher throughput of DNA sequencing and of general measurement systems using fluorescence-based imaging technology. This is a joint work with Bo Wang, Anqi Wang and Lei M Li.

An important test of a mathematical model is how well it describes a given phenomenon, which is often observed only indirectly and with error. This talk will provide an introduction to the statistical problem of inference for nonlinear dynamical systems. We will take a Bayesian perspective, constructing hierarchical models that incorporate prior knowledge and impose required constraints on model components.

The dynamics of two coexisting etiologic agents could be affected by factors such as cross-immunity, cross-enhancement, weather and vaccination. As an intent to analyze quantitatively the interplay of these factors, we describe a two pathogen model defined through a continuous time Markov jump process. We study the capacity of this model to sustain oscillations around its steady state and around its limit cycle as a function of these parameters. We use standard method to compute the power spectral density of the fluctuation of the stochastic model. Different scenarios are shown where vaccination may change the behavior of two pathogen. Our findings are applied to respiratory syncytial virus and influenza data.

In this talk I will give an overview of set-based methods for computing global dynamics. These methods combine topological methods, combinatorial and algebraic tools, and some rigorous numerics to build a collection of algorithms to compute a robust characterization of the global dynamics of a system. I will show examples and briefly describe software and some mathematical foundations.

Experimental data on gene regulation is mostly qualitative, where the only information available about pairwise interactions is the presence of either up-or down- regulation. Quantitative data is often subject to large uncertainty and is mostly in terms of fold differences. Given these realities, it is very diffcult to make reliable predictions using mathematical models. The current approach of choosing reasonable parameter values, a few initial conditions and then making predictions based on resulting solutions is severely subsampling both the parameter and phase space. This approach does not produce provable and reliable predictions. We present a new approach that uses continuous time Boolean networks as a platform for qualitative studies of gene regulation. We compute a Database for Dynamics, which rigorously approximates global dynamics over entire parameter space. The results obtained by this method provably capture the dynamics at a predetermined spatial scale. We apply our approach to study neighborhood of a given network in the space of networks. We start with a E2F-Rb network underlying the mammalian cell cycle restriction point and show that majority of the parameters support either the GO, NO-GO, or bistability between these two states. We then sample perturbations of this network and study robustness of this dynamics in the network space.

How many facial expressions of emotion can humans produce and visually recognize? How does our brain interpret them? Did some of these non-verbal signals evolved to serve as grammatical markers in human language? In this talk, I will overview the research we have done to address these fundamental questions. Specifically, we will show that a. people consistently produce at least 23 facial expressions of emotion (much more than the previously 6 proposed by Darwin), b. these emotions are perceived categorically, not as combinations of more basic categories, c. our visual system interprets them by analyzing shape and shading image features that represent the movement of individual facial muscles (known as action units, AUs), and d. three of these emotions have been compounded to form a (universal) grammatical marker of negation in human language. I will present a computational model of this perception of facial expressions and show behavioral and imaging experiments in favor of this mode. For example, we will identify a small region in posterior STS dedicated to the coding of individual AUs as predicted by the model. I will also show how this computational model can be used to define computer vision systems that yield superior results to state of the art algorithms and argue that the important problem of these vision tasks is detection, rather than recognition.

Contact structure impacts disease transmission on many different spatial scales (e.g. contact networks, mobility networks). Less well studied is how disease status impacts contact structure, and how this in turn impacts disease dynamics. In this talk I will look at this interplay from two different perspectives: first, in terms of a tradeoff between contact and transmissibility using a set of simple ODEs, and second, in terms of dynamics on a multi-layered network, where the structure of each later evolves in tandem with the disease. The first part is joint work with Chiu-Ju Lin and Kristen Deger, and the second part is joint with Karly Jacobsen, Mark Burch, and Grzegorz Rempala.

We consider a mutation-selection model of a population structured by the spatial variables and a trait variable which is the diffusion rate. Competition for resource is local in spatial variables, but nonlocal in the trait variable. We establish the existence and asymptotic profile of a steady state solution. Our result shows that in the limit of small mutation rate, the solution remains regular in the spatial variables and yet concentrates in the trait variable and forms a Dirac mass supported at the lowest diffusion rate. Furthermore, the steady state is unique and locally asymptotically stable. Global asymptotic stability is also established in a special case. Our result generalizes the finding of A. Hastings which says that in a two species competition in spatially heterogeneous environment, the slower diffuser always prevails, if all other things are held equal. This is joint work with Yuan Lou (Renmin Univ. and Ohio State).

Algebraic topology and dynamical systems are intimately related: the algebra may constrain or force the existence of certain dynamics. Morse homology is the prototypical theory grounded in this observation. In this talk we will discuss the computational Conley theory - a combinatorial-topological generalization of Morse theory, which is aimed toward capturing robust global dynamics across a wide range of parameters. The approach has by now been introduced multiple times during this emphasis semester, and our talk will focus on a capstone of the project - computation of the connection matrix. The connection matrix is the mathematical object which transforms the approach into a truly homological theory. We will discuss its significance and our current progress in its computation.

**This seminar will be held in room CH240 in Cockins Hall (1958 Neil Ave)**

A two-species biological depletion model in a bounded domain is investigated in which one species is a substrate and the other is an activator. Firstly, under the no-flux boundary condition, the asymptotic stability of constant steady-states is discussed. Secondly, by viewing the feed rate of the substrate as a parameter, the steady-state bifurcations from constant steady-states are analyzed both in one-dimensional kernel case and in two-dimensional kernel case. Finally, numerical simulations are presented to illustrate our theoretical results. The main tools adopted here include the stability theory, the bifurcation theory, the techniques of space decomposition and the implicit function theorem.

The network inference problem consists in reconstructing the structure or wiring diagram of a dynamic network from time-series data. Even though this problem has been studied in the past, there is no algorithm that guarantees perfect reconstruction of the structure of a dynamic network. In this talk I will present a framework and algorithm to solve the network inference problem for discrete-time networks that, given enough data, is guaranteed to reconstruct the structure with zero errors. The framework uses tools from algebraic geometry.

Collective cell movement is ubiquitous in biology, occurring both in normal processes (for example, development, wound healing) and in disease (for example, cancer). In most of these examples, how cells coordinate their movement is still not well understood. We will consider two examples: (i) angiogenesis is the process by which new blood vessels are formed in response to, for example, wounding or tumour growth. Typically, this has been modelled phenomenologically using the well-known snail-trail framework, leading to a coupled system of nonlinear partial differential equations for two key endothelial cell populations (tips and sprouts). Here, we revisit this model and show that a more formal derivation of the PDE model, from a discrete master equation framework, leads to a novel coupled system of PDEs to those studied in the literature; (ii) neural crest cell invasion is a very important early developmental process and also shares many common features with melanoma cell invasion. Here, we use a combined experimental and mathematical modelling study to shed light on a number of questions regarding the basic principles of this phenomenon.

The actin cytoskeleton is an important part of cellular motility. Expressions of actin regulators are altered in metastatic carcinomas. In this talk, I will discuss past and ongoing research on modeling two distinct motility structures involved in cancer metastasis, namely lamellipodia and invadopodia. Actin growth dynamics during the formation of the two structures show qualitatively similar trends, but closer analyses show distinct behaviors which suggests that key regulators may be playing different roles. In particular, one regulator called cofilin severs capped actin-filaments spurring further growth in the lamellipodia. Our mathematical model shows that its role during the formation of invadodia may be highly dependent on the availability of actin monomers.

Michales-Menten equation derived using the standard quasi steady-state approximation (sQSSA) has been widely used for the estimation of enzyme kinetic parameters. However, such approximation is only valid when enzyme concentrations is low. Thus, we found that the estimation can be biased when enzyme concentration is high. On the other hand, we found that a newly reduced model using the total QSSA (tQSSA) provides the accurate estimation of parameters regardless of enzyme concentration. For the estimation of parameters, we utilized the Markov Chain Monte Carlo (MCMC) based on the Bayesian approach to estimate parameter using non-informative and informative priors.

Symmetry methods have been used to great effect in solving problems in Mathematical Physics in the 20th century. With the rise of Mathematical Biology in the 21st century, it is only proper to consider the impact a group theoretic approach can make in this field. We will illustrate the usefulness of symmetry methods by considering two problems: a toy 2-strain infectious disease model and a chemotaxis model.

Due to the lack of kinetic information of many biochemical interactions, it is difficult to study and analyze continuous models of gene networks. The nonlinearity of these interactions makes the problem even less tractable. The Boolean framework (where gene activity is assumed to be either 0=OFF or 1=ON) provides a complementary approach that focuses on qualitative analysis. In this talk, I will show how the Boolean framework is used in modeling biological systems, and how tools from graph theory and combinatorics are well-suited for the analysis of these Boolean models. This talk will have the format of a tutorial.

Cell division is a highly dynamic stage of the cell cycle that requires both precise duplication of the cellular components as well as careful transport of the duplicated material to each daughter cell. Cells have to solve both a mechanical problem arising from the need to move replicated components to specific cell locations, and a timing and accuracy problem which is rooted in the need to make sure that the genetic material is correctly duplicated and that these processes progress within a well-defined time-span. These two problems are resolved in elegant ways across different cell types through the use of sophisticated nano-machines and biochemical reactions. We give an overview of mathematical modeling that describes mechanisms by which cells not only move components but also control their positioning and timing so as to generate well defined protein patters that can then be used as spatial cues for division. We discuss simple mechanisms for intra cellular pattern formation and apply them to FtZ ring placement mechanisms in some dividing bacterial cells.

Networks which show the relationships within and between complex systems are key tools in a variety of current scientific areas. A central aim in network analysis is to find a suitable metric for network similarity and comparison. I'll discuss a definition for the space of all networks, and show that this definition leads to natural notions of (1) isomorphism between pairs of networks, and (2) a compatible distance between networks. Then I will discuss the computational complexity involved in computing the network distance, and develop lower bounds by using invariants of networks that are significantly simpler to compute. Finally, I will describe a family of invariants which we call "motif sets" that characterize networks up to isomorphism.

The problem of understanding stable steady state patterns appears in processes such as cell differentiation, memory formation, and multistationarity of gene networks. Systems of differential equations have been used to study these problems, but the size and nonlinearity of the systems make the problem difficult. In this preliminary report we will present a combinatorial approach to study the steady state patterns of differential equations using AND gates. Our results show that if the interactions between variables in a gene network are made up of AND gates, then we can find the stable steady state patterns of a network of arbitrary connectivity and size.

Stochasticity arising from low copy numbers of chemicals in reacting systems can have a huge impact on the results of numerical simulations. However, such strong noise is often very difficult to analyse. As the population of the reactants grows we expect the influence of randomness to decrease. This puts us in the "weak noise" regime, where noise still has an effect, but the dominant contribution is from the mean-field kinetics. In this talk I review the techniques that can be used to understand this weak-noise regime from a spatio-temporal patterning point of view.

Being able to create and sustain robust, spatial-temporal inhomogeneity is an important concept in developmental biology. Generally, the mathematical treatments of these biological systems have used continuum hypotheses of the reacting populations, which ignores any sources of intrinsic stochastic effects. We address this concern by developing analytical Fourier methods which allow us to probe the probabilistic framework. Further, a novel description of domain growth is produced, which is able to rigorously link the mean-field and stochastic descriptions. Finally, through combining all of these ideas, it is shown that the description of diffusion on a growing domain is non-unique and, due to these distinct descriptions, diffusion is able to support patterning without the addition of further kinetics.

The extracellular matrix (ECM) provides rich biochemical and biomechanical cues that regulate many cell processes such as adhesion, migration, proliferation, differentiation, apoptosis, and intercellular signaling. Cells continuously maintain their local ECM through a balance of matrix synthesis, chemical modification, degradation, and the application of traction forces. It is believed that a breakdown of this balance and subsequent changes to the mechanical properties of the ECM could disrupt homeostasis and promote pathogenesis such as cancer metastasis and fibrosis. The dynamic mechanical interactions of cells with ECM are modeled in a physically realistic agent-based model. Resulting from these simple mechanical rules, simulated cell-populated ECM exhibit many emergent behaviors observe in experimental systems including matrix remodeling and directional cell migration.

The study of the development of the embryonic neural tube under the patterning action of various morphogens constitutes an important research topic within Developmental Biology. The use of multiscale models based on a combination of ordinary and partial differential equations is a well established research paradigm in this area by now. Linear diffusion equations are widely used to describe morphogen propagation. However, this modeling assumption is challenged by a number of experimental results. As a proof of concept, we try to circumvent these difficulties by introducing and analyzing a nonlinear diffusion model yielding finite speed of morphogen propagation in the context of dorso-ventral patterning of the neural tube.

This talk will review various stochastic models and algorithms for reaction-diffusion systems with multiscale nature, which have been previously developed. Then, I will introduce our algorithms for stochastic reaction-diffusion systems coupling different modeling schemes for better efficiency of the simulation. We can apply the algorithms to the systems including regions with a few molecules and regions with a large number of molecules. A domain of interest is divided into two subsets where continuous-time Markov chain models and stochastic partial differential equations (SPDEs) can be respectively used. The talk will present simulation results of several examples. This is a joint work with Radek Erban at the University of Oxford.

The expression of a gene is usually controlled by the regulatory elements in its promoter region. However, it has long been hypothesized that, in complex genomes such as that of the human, a gene may be controlled by distal enhancers and repressors, not merely by regulatory elements in its promoter. Globally, the chromatin structure is closely linked to the biological function of the genome, and has great implications in human health and disease. A recently developed molecular technique is able to detect physical contacts between distant genomic loci, validating the theory that communications between such elements are achieved through spatial organization (looping) of chromosomes to bring genes and their regulatory elements into close proximity. Such a molecular technique, coupled with the Next Generation Sequencing (NGS) technology, enables genome-wide detection of physical contacts. The availability of such data makes it possible to reconstruct the underlying three-dimensional (3D) spatial chromatin structure of a genome and to study spatial gene regulation. However, several special features of the NGS-based high-throughput data have posed challenges for statistical modeling and inference. In this talk, I will focus on discussing strategies for addressing such features, including dependency, sparsity, and over-dispersion. Examples will be provided to illustrate the importance of taking such features into consideration.

Malignant gliomas are the most common type of brain cancer, which arise from glial cells, and in their most aggressive form are called GBMs. GBMs are highly invasive and difficult to treat because cells migrate into surrounding healthy brain tissue rapidly, and thus these tumors are difficult to completely remove surgically. GIMs, which can comprise up to one third of the total tumor mass, are present in both intact glioma tissue and necrotic areas. They apparently originate from both resident brain macrophages (microglia) and newly recruited monocyte-derived macrophages from the circulation. Activated GIMs exhibit several phenotypes: one called M1 for classically activated, tumor suppressive, and another called M2 for alternatively activated, tumor promoting, and immunosuppressive. Within a tumor the balance between these phenotypes is typically shifted to the M2 form. Numerous factors secreted by glioma cells can influence GIM recruitment and phenotypic switching, including growth factors, chemokines, cytokines and matrix proteins. In this work, we focus on mutual interaction between a glioma and M1/M2 microglia mediated by CSF-1, TGFbeta, and EGF. Up-regulated TGFbeta leads to up-regulation of Smad within the tumor cells and secretion of MMPs, leading to proteolysis for EMT process and cell infiltration. The mathematical model consists of spatial and temporal evolution of densities of glioma cells, M1 type cells, M2 type cells, and concentrations of CSF-1, EGF, TGFbeta, Extracellular matrix, and MMPs in the framework of partial differential equations. We developed the model to investigate the mutual interactions between tumor cells in the upper chamber and microglia in the lower chamber. In the experiments, Boyden invasion assay was used to show that this mutual interaction induces glioma infiltration in vitro and in vivo. We show that our simulation results are in good agreement with the experimental data and we generate several hypotheses that should be tested in future experiments in vivo. This interaction between a glioma and its microenvironment suggest a great potential in oncolytic virus (OV) therapy because of its complexity and immune suppression or responses. We discuss recent development in OV therapy for better understanding of the whole system and suggest new strategies to eradicate glioma cells in the framework of various combined therapies.

**This talk will be held in Jennings Hall, room 360, and will not be streamed online.**

Nematode worms of the genus Caenorhabditis rely on an egg laying structure, called the vulva, to properly expel embryos from the adult body. Development of the vulva involves six vulval precursor cells (VPCs) adopting one of three possible cell fates: primary, secondary and tertiary. The signaling pathways involved in vulval development are highly conserved across nematode worms, although data from our group shows differential responses to partial pathway disruption between two species, C. elegans and C. briggsae. In order to investigate the pathway components that may be responsible for these species-specific differences, we have developed a biologically based mathematical model of the main signaling networks (EGF/Ras, Notch and Wnt) responsible for VPC cell fate specification. We have identified key components in the model that account for the observed species-specific responses. In particular, we find that C. elegans-like parameter sets are more sensitive to disruptions in the EGF/Ras signaling pathway than C. briggsae-like parameter sets, and that differences in the spatial distribution of Wnt signaling may be responsible for the C. briggsae robustness. Taken together, the theoretical results have allowed us to identify modifications in the underlying interaction network that may give rise to the experimentally observed phenotypic variability.

TBD

We considered a combination therapy of cancer. One drug is a vaccine which activates dendritic cells so that they induce more T cells to infiltrate the tumor. The other drug is a checkpoint inhibitor, which enables the T cells to remain active against the cancer cells. To address the question of synergy between the two drugs, we develop a mathematical model using a system of partial differential equations. The variables include dendritic and cancer cells, CD4+ and CD8+T cells, IL-12 and IL-2, GM-CSF produced by the vaccine, and a T cell checkpoint inhibitor associated with PD-1. We use the model to explore the efficacy of the two drugs, separately and in combination, and compare the simulations with data from mouse experiments. We next introduce the concept of synergy between the drugs and develop a synergy map which suggests the proportion in which to administer the drugs in order to achieve the maximum efficacy.