Dynamical and Mechanistic Reconstructive Approaches of T Lymphocyte Dynamics: Using Visual Modeling Languages to Bridge the Gap between Immunologists, Theoreticians, and Programmers.

Dynamic modeling of lymphocyte behavior has primarily been based on populations based differential equations or on cellular agents moving in space and interacting each other. The final steps of this modeling effort are expressed in a code written in a programing language. On account of the complete lack of standardization of the different steps to proceed, we have to deplore poor communication and sharing between experimentalists, theoreticians and programmers. The adoption of diagrammatic visual computer language should however greatly help the immunologists to better communicate, to more easily identify the models similarities and facilitate the reuse and extension of existing software models. Since immunologists often conceptualize the dynamical evolution of immune systems in terms of "state-transitions" of biological objects, we promote the use of unified modeling language (UML) state-transition diagram. To demonstrate the feasibility of this approach, we present a UML refactoring of two published models on thymocyte differentiation. Originally built with different modeling strategies, a mathematical ordinary differential equation-based model and a cellular automata model, the two models are now in the same visual formalism and can be compared.

The perspective is to encourage immunologists involved into mathematical modeling or software productions, to adopt a visual graphical language, here mainly the unified modeling language (UML) “state-transition” diagram to ease the communication, the reuse and the extension of their models.

Complexity of the Immune System

The immune system is a complex biological adaptive, highly diversified, robust and resilient system, characterized by complexity at different levels. Lymphocytes are the central actors of the immune system, in the middle of a multi-scale biological organization, “from molecule to organism”. Multi-scale modeling remains a challenge, as for other biological systems (1). Despite recent systems biology initiatives to understand and model the immune system (2), we are still far from having the appropriate tools to understand its dynamics and to easily communicate among various researchers who observe this system at different levels of granularity and attempt through software modeling to answer different questions. Several complementary experimental methods and models have been used to explore lymphocyte dynamics and turnover (3, 4) and to model it in health, aging and diseases (5).

Drawbacks of Current Dynamics Lymphocyte Modeling and Evolution

System dynamics models deal with time, formalized with two distinct concepts as “discrete time,” by a succession of time points and intervals, or as “continuous time” (6). Models of lymphocyte population dynamics and turnover (3) have primarily been based on mechanistic reconstruction with continuous time models. The fluxes of cell populations are then described by differential equations. These mathematical models describe for example the thymocyte differentiation and selection (7), until the thymic export (8–10), the homeostasis of CD4/CD8 T cells (11, 12), the CD8 immune response (13, 14), or the Bromodeoxyuridine or deuterium labeling (15) to account for turnover. Up to now, simulations and validation of some of these models reveal interesting T cell dynamics properties: how the system grows, self-maintains as well as the effects of perturbations, i.e. how the system reacts to antigens, collapses and reorganizes. However, integrating the heterogeneity of cell populations, phenotypes, lineages, cell location and interactions, cell differentiation across generations (16) in the different biological, and time scales, is problematic in such a mathematical form, which make these models particularly difficult to handle.

The evolution of homogeneous mathematical model of cell populations (17) toward “spatialized,” discrete, and heterogeneous software models (18) has allowed the reproduction and observation of more detailed and thus complex behaviors. For example, this made possible to model lymphocyte dynamics from thymic selection (19, 20) up to quantitative modeling of immune responses, as extensively reviewed (21) with development of agent-based and automata models (22). However, both population-based mathematical model (a top-down approach) and discrete cell-based model (a bottom-up approach) and the various platforms developed have limitations (23). Conversely, the Statecharts language (24) and the visual reactive tools (25) such as biocharts (26) and reactive animation applied to various systems (27) developed by Harel et al. are a powerful way to simulate complex dynamical biological behavior with more didactic representation than equations. Such models have revealed emergent properties during thymic differentiation (28) pancreatic islet organogenesis but also the immune response in lymph node (29).

Lack of Interoperability, Under-Use of Software Models

Although in Immunology there is more than 20-years tradition of software and mathematical modeling, very few of them have been the object of further exploitation once published and made available (30). Models are often under-used because experimentalists can be reluctant to entertain mathematical formalization and because published models are largely disposable: rapidly forgotten after being published, instead of providing a foundation to build upon. Moreover, the various expressions of these models with different mathematical descriptions, programing languages, software libraries and graphical packages, require much effort in understanding and running the software and prevent interoperability.

Using Visual Language to Communicate and Execute Models

Immunologists often conceptualize the dynamical evolution of their systems in terms of “state-transitions” of biological objects and do it by means of personalized and informal graphical illustrations. Thus, the adoption of a more formal and standard type of state-transition diagram could improve the current situation to not only help biologists to better understand each other but also to facilitate the production and the reading of software code executing these visual transitions, at level of populations or agents (31).

Thus, in this paper, we promote the development and usage of a visual, computational approach more comprehensible than mathematical equations and programing instructions. This should improve our understanding of lymphocyte dynamics, the exchange on this understanding and simplify the implementation of models by non-specialists delivered from the production of executable code or mathematical equations, to concentrate to in silico experiments.

Level of abstraction and multi-scale modeling

A model describes a complex system from the “real word” and thus requires abstraction. This abstraction is performed as the immunologist decides about an experimental protocol in order to observe selected objects at different scales and to follow them in time and space. For example, the capacities of the immune system to preserve the homeostasis and to provide rapid adaptation to an antigen and anamnestic responses can be observed at the organism level, through physiological or pathological clinical observations that relate to lower scale levels. At molecular level, the somatic generation of the diversity of an immuno-receptor, as the TCR, allows for a dynamic network of interactions with antigens. At the cell level, this leads to clonal selection, activation and division. At the organ level, the fluidity of the system insures constant tissue redistribution of cells and molecules, cell migration from thymus to spleen and lymph nodes.

Thus, models of lymphocyte population dynamics and turnover consist in reconstructing the components or “entities” of the system across various scales, from molecules to organisms, to determine the relations/interactions through space (varying from micrometer to meters) and “processes” through time (varying between microseconds to years) as explained below. However, the formalization and abstraction of the immune objects as entities undergoing processes, with the help of spatial and dynamic ontologies, respectively defined as SNAP and SPAN (32), as well as cell/molecule interactions (33), is rarely done, maintaining a language-barrier between biologists and theoreticians. In the following, some examples will be given to help the immunologists with the transition between current mathematical models to computer ones and with the terms currently used in modeling.

Define entities, states, location, interactions, granularity

The immune “entities” could be described according to the language used by the modeler. A cell exists in one “state”: it could be quiescent or in a given phase of the cell cycle or dead. In addition the phenotype and/or a function of a cell define a given state, as CD4 helper T cells. Cells are “located” in various tissues and are in “relation” with other entities. Finally, cells can be considered at various level of “granularity.” For example, T lymphocyte populations are “aggregation” of T cells according to criteria of phenotype, structure or function, although heterogeneity still prevails inside these populations at lower granularity. Accordingly, cells can be modeled at population level (with continuous model as ordinary equation) or at cell level according to space (with discrete model as automata or multi-agent system).

Define processes

According to ontologies, cells participate to various processes, such as division, activation, differentiation, interaction, clonotype selection, apoptosis or migration. According to the states of the cells, their evolution can thus be modeled as “state-transition” that can be applied to various processes in parallel: for example, a thymocyte can differentiate while migrating in cortex and medulla. Note that processes at other levels like molecular or organ levels can similarly be described and modeled. Finally, all these process will determine the global cell dynamics and turnover.

The Unified Modeling Language for High-Level Modeling

“High-level programing languages” are based on abstraction and use of natural language that is easier to understand as compared to “low level programing language,” based on codes. Thus, visual modeling language considers biological-object as conceptual abstract-objects that endure processes. The level of abstraction allowed by these diagrams makes possible to distinguish more easily the “entities” as T cells and the “processes” that occur at different levels such as differentiation, migration and cell cycle. Moreover, such “state-transition diagrams” allow computing parallel pathways at various scales to avoid redundancy that is inherent in the formal description of multi-level, heterogeneous and concurrent systems and to model heterogeneity in a very simplified and economical form (as compared to mathematical equations) (31).We have thus used the well-established Unified Modelling Language (UML – a sofware standard) that still remains approachable to the lab-immunologist, convenient for the theoretician and that can be directly adopted for the high-level graphical depiction (31, 34). The adoption of UML state-transition diagrams that transcends any programing language or computer platform, will allow both experimentalists and theorists to work together at a higher level than writing software code or mathematical equations. This final step is progressively more and more automatized out of the diagrams. Example of basic transformations of mathematical equations into state-transition diagrams including elementary, parallel, independent, or coupled state transition have been given (31), to familiarize biologists with the general approach.

Refactoring from Low Level (Code, Equations) to High-Level (Diagram) Modeling Languages

To convince the immunologist of the feasibility of this approach as well as the benefit gained by adopting it, we sketch in the rest of the paper how existing “low level programed” immune models should gain in readability and accessibility by adopting a “high-level graphical” representation under the form of “state-transition diagrams.” We present a “refactoring” of two published models of T cell biology in the thymus. Refactoring consists in restructuring the code or equations of a model to improve its expression, readability and extensibility, without changing its external behavior. One model consists in cell population differentiation modeling with differential equations (continuous model). The other one is a discrete model. Originally it was an automata model consisting of a discrete lattice, where each site (cell) in a given state, follows some rules in space and time that depends on local neighbors (18). It has been refactorized as an agent-based model (ABM), depicting individual cell behavior through thymus differentiation and migration. It would be much too long and redundant to describe in details the behavior of these two models. We do not pretend here to modify at all the results obtained by the running of these models (the readers interested in these results are invited to access the original papers). We have just reshape them into a state-transition diagrammatic form that allows execution of simulations reproducing the original results with similar parameter values.

Population-based model describing the conveyor-belt T cell differentiation in thymus

The original model (8) is a compartmentalized ordinary differential equation (ODE) model, rather complex to read and manage by immunologists. This model reflects the conceptual “conveyor belt” model of thymic T cell differentiation, schematically represented by immunologists by the continuous ordered transition of cells through the different stages with time (35–37). Figure 1 represents a biological schema, originally published and the “state-transition” description of the model (in a UML state-transition diagram) as it is proposed now. Although the original model is composed of 30 differential equations, the whole mathematical description and the code that captures it, can easily be deduced and regenerated from the Figure 1. Conversely, the mathematical equations can be automatically generated from the state-transition diagram as previously described (38). In essence, the model is summarized by the input, transition and output from the thymus, by “parallel processes” that occur concomitantly, as differentiation, cell cycle, proliferation and death, and by exit from the thymus. Note that these parallel processes concern various biological levels and time scales. The “differentiation” process represents each stage of the conveyor belt, from double negative (DN) to single positive (SP) cells, with flows into and from a particular stage according to the general equation: p, d, and u represent proliferation, death, and differentiation, respectively, and x represents the ith stage. The model mainly consists of constant hematopoietic progenitor influx in thymus (Sn); differentiation between thymocyte developmental phenotypes DN and double-positive (DP) cells differentiation into CD4+ or CD8+ cells, then egression of SP stage, either, to the periphery [Us4(i), Us8(i)]; proliferation [Pn(y), Pp(y), and Ps(y)]; positive and negative selection (a4, a8); and natural cell death (Dn, Dp, and Ds). In parallel, the “cell cycle” is represented: the cell switches between quiescence (G0) and cycle with division (S/M). The parameter γ set to 1 represents the cell division into two daughter cells. There is the possibility to induce a perturbation into the system through the specific depletion of T cells entering the S/M cycle phase (p), if γ is set to 0. This represents the presence or absence of a pharmaco-genetic conditional treatment by ganciclovir that induces apoptosis related to the incorporation of a nucleotide analog during DNA elongation. This rule applies to all cell populations except in late DP quiescent cells as indicated in the schema. The “proliferation” depicts that the daughters of a proliferating cell transit into the next generational compartment, except during treatment (γ = 0) when dividing cells die by apoptosis and are lost from the model. The parameter u is an increasing function of generation (G1 to n), making cells more likely to differentiate between phenotypic compartments as they progress through the cell generations.

Figure 1

Refactoring the continuous population (ODE) conveyor-belt model of thymocyte differentiation (. This diagram represents the evolution of cell populations in the thymus, represented by their state from DN to SP, with UML figuration of the input (close circle arrow) and output representing death and thymus exit (open circle arrow) and transitions (oriented arrows). Parallel processes (underlined with red box) as differentiation, cell cycle and proliferation are depicted more explicitly in this representation than with the original 30 mathematical equations. Annotation with the proliferation and death rates values indicated in each state are values from the best scenario observed in the original paper.

The parallelism in this graphical model largely simplifies the original formulation of ODEs while remaining faithful to it. Hierarchy and compound states are present again clearly reducing the diagram clutter. Other representation of the model depicting the differentiation with linear cell generations is also possible (31).

Discrete model describing the differentiation along the migration of T cells in the thymus in a 2-D environment

The original model (20) is a discrete-based “cellular automata” computer model. The model depicts the behavior of individual thymocytes that evolve in the 2-D epithelial cell network, guided by chemokines gradients. The current model (Figure 2) is now an Agent-Based Model (ABM). Again, the interested reader is referred to the original paper for a detailed understanding of the simulation. Although available for download, the 40 pages of FORTRAN source code are far from easy to understand. After refactoring, the transition rules of any agent (thymocytes) map onto a parallel state-transition diagram. The conception of the state-transition diagram, as done here, should considerably improve the understanding of the model (even for the original programmer) allowing the researchers to progress further with the existing simulator. A complete description of the model includes additional implementation details that are abstractions of the mechanisms behind cell decisions to differentiate. The parallel state-transition diagram represents the different simultaneous transitions taking place in the model and coded as various cellular automata rules: a cell in the model as it differentiates transits in successive states, it may be bound or not to thymic epithelial cells via TCR/MHC, it moves and may be located into one of several anatomical compartments of the thymus. As indicated in the boxes, the gradient of chemokines (k) orients the migration of cells for each specific stage of differentiation, and chemokines are localized in specific areas of the thymus. A T cell sums up the time and number of interactions with the same or different epithelial cells. This sum value determines whether the T cell is positively or negatively selected. DP and SP phenotypes have their own threshold parameters. They cannot differentiate until they cross a threshold number of interactions. If, after a given time, a DP cell has not reached this threshold, it enters apoptosis by neglect. If the time is too long, it is negatively selected. A threshold parameter simulates the phenotype decision. With this “signal-duration” hypothesis, long duration TCR-MHC interactions promote the CD4 phenotype and short duration promotes the CD8 phenotype.

Figure 2

Refactoring the automata single-cell model of thymocyte differentiation and migration (. The 40 pages of Fortran code is transformed in computer UML state-transition diagram. This diagrams describes the experimentally observable heterogeneity, and the biologically relevant parallel processes (underlined in red box). As in Figure 1, input, output, and oriented transitions are described in the state-transition diagram. Parallel state-transitions represent the evolution of single cells in the thymus, represented by their differentiation from DN to SP stage, sequential binding event of TCR/MHC peptide on epithelial cells, thymic location, egression of cell when matured. Additional qualitative abstractions for computational model of individual cells in semi-realistic environment are represented by the 2-D array: cells are allowed migrating sequentially through the epithelial cell network (black network) across the various thymic areas, guided by chemokine gradients CXCL12 (red), CCL19/CCL21 (green), S1P (blue).

It is important to notice that both refactored models, population-based and ABM can now be compared, are directly executable and can provide simulations of physiology, pathologies, and treatment, while not being the scope of this paper. Moreover, their parameters can be automatically tuned to fit experimental data. Any ABM model can also be run as a population version to save time in simulation (McEwan et al., manuscript in preparation). Moreover, the flexibility of these diagrams allows assembling the parts of the biological puzzle piece after piece and improving the models.

Perspectives

As shown here, state-transition diagrams can represent high-level semantics suitable to clarify immunological concepts and to aid communication among interdisciplinary researchers. It can also represent low levels quantitative information suitable for individual-based ABM and population-based ODE modeling. Organization of immune knowledge using a standardized, diagrammatic formal language should greatly improve knowledge integration at multi-scale levels and sharing between experimentalist and theoretician collaborators, rendering their software more readable, scalable, and usable. We are currently working on ways of automatically generating executable code out of these state-transition diagrams. State-transition diagrams supports the extension and interoperability of published models. This will help for dynamic computational modeling of lymphocyte behavior in health and diseases, and for “in silico” experiments to predict and explain the puzzling T cells dynamics and the effect of immunological perturbations.

Author Contributions

Véronique Thomas-Vaslin and Hugues Bersini designed the project and wrote the article, Adrien Six and Jean-Gabriel Ganascia participated to discussions.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.