Multiscale Modeling and Simulation of  The Immune System

Shlomo Ta'asan
Carnegie Mellon University



Existing mathematical models

The complexity and dynamic nature of the immune system stimulated many mathematician/modelers that attempted to gain understanding into its functioning and regulation. Indeed several mathematical models of different aspects of the immune system have been developed. The translation of the immune network theory of Niels Jerne (64) into mathematical models has been one of the first important developments in the mathematical modeling of immunology (65-69). Immunologists enthusiasm to this line of research probably faded away at the time the clonal selection theory was introduced. Mathematical models have proved very useful in the study of some aspects of the dynamics of AIDS (7,8), in particular, in interpreting the effect of treatments (10,11,23,24), as well as in certain aspect of the disease and its progression (12-22,25-27,29), and some aspects which have not yet been understood, for example, the latent phase (9,28). The study of T-cell activation and of the cognate interaction (46,47,48,50,79) has benefited from the use of mathematical models (44,52,53)). They are more quantitative and as a result allow for a more precise and refined analysis of the dynamics of interactions of the receptors and the activation of the cells of the immune system (55,58,59). More recently simple mathematical models were used to discuss the immune memory, Anita (91).

None of the aforementioned approaches were designed to model the biological aspect of the immune response. Computer based approaches like cellular automata or complex system are the only approaches existing so far, which have the potential to describe adequately such processes. In the context of the immune system, there have been at least two lines of research following that approach. P. Seiden and F. Celada have used cellular automata to capture such processes, (86,87,88,89). The goal was to capture as much as possible of the dynamics of the immune system and to be able to make experiments ex machina. A study of the thymus along similar lines has been developed (77,78). This interesting approach has provided important insights on the regulation of the positive and negative selection, and on the dynamics of production of the TCR repertoire in the thymus.

Continuous models

Many phenomena exist on a wide range of scales and the proper modeling tool is sometimes dictated by the question of interest. While it is always possible to use a detailed microscopic model, is not an economical approach for dealing with macroscopic dynamics. A macroscopic model in such a case, which usually use continuum approaches is best appropriate. Different problems in physics and chemistry have been successfully treated using continuum approaches, i.e. either ODE or PDE. However, numerous examples exist in which the phenomena is understood much better on a microscopic level, and a good insight into the dynamics at larger scales is missing. For example, the behavior of polycrystalline materials which exhibit complex behavior due to microstructures in them; a lot is known about the dynamics at the micro level, while quantitative models for the macroscopic behavior are far from being satisfactory. Construction of macroscopic models directly in such cases is not feasible, and the micro-level behavior has to be integrated in some form into the macroscopic models.

Partial differential equations as well as ordinary differential equation can be viewed as limits of discrete models. The continuum equations of fluid dynamics or of elasticity can be explained as limiting behavior of atomistic/molecular models as the number of molecules/atoms is large enough. These limiting processes are valid under certain smoothness assumption for the spatial distributions and the presence of enough 'particles'.

Many biological systems involve high levels of heterogeneity and variability, features that raise questions regarding the appropriateness of continuum models. Consider a classical reaction equation between a site S and a ligand L  It is a classical result that the concentration, say S, obey an ordinary differential equation. This equation describing the dynamics of concentration is valid under the assumption that the mixture is homogeneous in space, that enough sites S and ligands L exist, and that there is no variability in each type. Considering biochemical reactions happening as part of the immune response, these assumptions are definitely incorrect. The concentration of antigen, for example, is not uniform, and similarly for almost any other population relevant to the immune response. This suggests the use of partial differential equations (PDE) instead of ordinary differential equation (ODE), and in some case it is justified. However, when spatial distributions are not smooth enough or at low levels of concentration even a PDE model breaks down. These observations are relevant when dealing with the immune system. The B-cell population relevant to a given antigen, their variability in terms of affinity to the antigen certainly do not satisfy the requirements for a continuum description, at least at the initial response. Moreover, during the immune response a subpopulation of B cells goes through a significant mutation phase and a single spatial variable cannot capture the right dynamics.

Discrete models and their continuous limits

Certain aspects of biological systems are better captured in a discrete models which allow for more generality than continuous models even under the restriction of the dynamics to lattices.  The main feature in biological systems which is different from most chemical or physical systems is the high level of heterogeneity and variability, which makes the statistics different. Although we deal with large numbers, the large variability in the elements 'state' require statistic on not-so-large numbers. Simulation on the other hand have to be involved with large numbers to reflect the variability in the system. Continuum model which are common in other sciences may need to be augmented by stochastic 'particle-like' models when dealing with biological systems. In certain cases, but not all, the limiting behavior reduces to continuum models.

Coarse-graining and local equilibrium

The choice of modeling approach is closely related to coarse graining questions. The route to macroscopic dynamics goes through distribution functions. Most large scale interacting system have a common feature allowing the passage to macroscopic description, i.e., the existence of local equilibrium functions. These are distribution functions describing some aspect of the problem, and have the tendency to approach their limiting behavior quite fast regardless of initial distribution. This is the result of the large number of interactions in the system, which is responsible for the relatively simple dynamics of the larger scales (a dynamics that can be described by a few macroscopic variables). The precise form of these distribution functions determine the form of the macroscopic equations, and therefore the dynamics itself. For example, molecules in a gas admit a velocity distribution that is bell shaped (Maxwellian), and this is obtained after a relatively short time independent of the initial velocity distribution. A purely Maxwellian distribution leads to the Euler equation of fluid dynamics, while a small perturbations from this distribution leads to the Navier-Stokes equations. Thus, the study of these functions may lead to insight into the emergence of the macroscopic picture.

In terms of the immune system such functions may not be not known but they presumably exist. This may be relevant, for example, to distribution functions for the level of expression of cell receptor population. These distribution functions have a single attractor in simple systems. In immune system, due to its complex response patterns, it is likely that a few attractors exist. The evolution of these distribution functions plays a major role in the global dynamics. A possible scenario is that the wandering of the distribution functions among the limit points can be an explanation of the change in character (Th1/Th2), bifurcation or chaos. Therefore, part of the more theoretical study in this work will be geared toward understanding these attracting limits and their role in construction of a macroscopic model.

Experiments and modeling

The immune response to pathogens is a complex, highly regulated system involving numerous interactions between different cell types. The cells of the immune system communicate with each other by direct cell-cell contact and deliver signals to each other directly, through cell surface molecules, or indirectly, via secreted proteins, known as cytokines. One of the challenges is to understand the integrated behavior of the innate and the adaptive systems. It is well agreed that the two work in concert; the adaptive system is triggered by the innate response, it is modulated by it, and ultimately it uses members of the innate system to eliminate pathogens. Separation of the innate and adaptive system is therefore undesirable if we want to understand the full system behavior.

The tight control over the recognition of foreign pathogens and development of appropriate effector functions is poorly understood at the present time. This is partly due to the fact that, in order to understand how and why a particular response develops, it is becoming necessary to devise experiments that can look at the whole immune system rather than its individual parts. In the past, experimental techniques that isolate the questions of interest have been used. For example, the actions of a particular cytokine can be determined using in vitro experiments with individual cell types and purified cytokines. In addition, the in vivo relevance of the same cytokine can be determined by generating a mouse lacking this cytokine through the targeted disruption of the cytokine gene. These approaches have been extremely successful in determining the functions of individual cytokines or cell surface molecules. These experiments do not, however, allow a true assessment of the relative importance of each molecule in the naturally occurring immune response, in which all of the cells and cytokines are functioning simultaneously. A good example of this is the case of interleukin (IL)-2 which, from in vitro experiments, was thought to be the most important factor driving T cell proliferation (1, 2). It was a surprise, therefore, that IL-2 knockout mice showed little to no impairment of T cell proliferation, but rather appeared to have defects in T cell death, which led to the development of inflammatory bowel disease (3-5). Do these results mean that IL-2 is not an important factor for T cell proliferation in vivo? Probably not: it is likely that, in the normal mouse, IL-2 does play an important role in T cell proliferation, but that the expansion of pathogen-specific T cells is such an important function of the immune system that other mechanisms can be utilized in the event of IL-2 dysfunction.

These, and other similar experiments, point to the flaws in using reductionism experimental techniques to address complex issues of immune regulation and development. There is now an increasing need to understand how the immune system functions as a whole in order to be able to predict outcomes of a particular challenge better, and to potentially devise strategies that could influence particular immune responses. An important concept is that the immune response evolves rapidly in time and its interactions are all highly regulated.

A major difficulty in gaining an insight into different processes in the immune system has to do with the limitation on possible measurements that can be conducted. Measurements affect the system and may hurt its integrity. The most accessible measurements are those involved with blood samples. Although these are easily done, the number of blood samples from a given mouse is limited as it has only 2 ml of blood. Thus, measurements of dynamic processes are very limited. In addition, there are many factors that cannot be measured without killing the animal involved; detailed examination of lymph node which is essential to understand the adaptive response. These and other limitations on the availability of experimental data regarding the system pose challenges to the immunologist that can be overcome only by the use of modeling together with system identification techniques. This is one of the places where the scientific questions of immunology origin are tightly connected to mathematical questions. Given a set of measurements of a certain type, what can we say about the internal structure of the system? How stable is the identifiability process? What is the relation between errors in measurements and errors in parameter estimation? Such questions have been studied in the context of other problems in various sciences. There are mathematical frameworks to deal with them as well as numerical algorithms for fast calculations. We regard this part of the work as a very fundamental to the development. It is only by a tight connection between experimental data, modeling and system identification that one can get hold of the non-measurable aspects of the system, and enables a deeper understanding of the immune system.

There has been relatively little interest in this area, on the part of experimental immunologists, probably due to the complexity of the mathematics necessary for this work and the dearth of mathematicians sufficiently cognizant of the immune system to make such models useful.

Our Modeling Approach: Continuous-discrete multiscale approach

Our proposed direction of research is guided by the scientific issues mentioned before, and the approach proposed here answers all the major difficulties presented. The most closely related previous work was done by Seiden and Celada (87-89), and parts of our model resembles their work. However, the modeling approach presented here is different in several respects, as we mentioned before: 1) the model involves {\em all} the subsystems; 2) it uses multiscale representation; 3) it involves a hybrid discrete-continuous representation; 4) it invokes parameter estimation methodology for quantitative prediction capabilities.

Among the unique features of this work is the inclusion of both innate and adaptive systems. This by itself is not enough to tackle the complexity of the system. Mathematical tools which enable us to make connection between models of different types (discrete vs. continuous) will be used; stochastic model, ODE's and PDE's which have been used successfully for different purposes are related to each other through certain limiting processes. Understanding their relation is probably more crucial in biological problems than in other fields due to the natural heterogeneity in such systems.

Our model is stochastic and include the other modeling approaches, as well as specific models using ODE/PDE as limiting cases where certain parameter (population size, for example) is sufficiently large, and certain uniformity exists. Moreover, our model tackles the problem of multitude of scale present in the problem and suggests an effective way to handle it, thus, it is not a mere extension of previous work. While previous work focused on important issues such as affinity maturation, clonal expansion etc, the relation to real experimental data was not emphasized enough. Different models can give qualitative pictures that resemble the experimental data, however, without a quantitative comparison one cannot rule out the possibility of missing crucial factors in the model. Our effort will therefore try to bridge this gap using appropriate techniques described later.

 A general way to accelerate computationally intensive simulation is to identify different scales in the problem (spatial and temporal), leading to multilevel approaches. This has been applied to many problems of different types, including solutions of PDE's, optimization techniques, integral equations, particle simulations and more. Multiscale techniques answer both the problem of efficient representation as well as the fast algorithm issues. This was a guideline for us in choosing the proposed approach. The main model is that of a mesoscale level but its precise definition comes from interactions with subsystem models. This is needed in order to have an accurate representation of the complex interaction that take place in the immune system. It allows the tuning of certain parts of the mesoscale model as well as to have the right dependencies between different parts whose behavior depends on some smaller scale 'entities'. A typical connectivity of such models is given in figure 1 where three levels of modeling are shown. We believe that at least two will be needed for accurate model, although a single level on the mesoscale will already give new information.