# Modeling Early Stages of Bone and Joint Infections Dynamics in Humans: A Multi-Agent, Multi-System Based Model

Salma Alsassa, Thomas Lefèvre, Vincent Laugier, Eric Stindel, Séverine Ansart

► To cite this version:

Salma Alsassa, Thomas Lefèvre, Vincent Laugier, Eric Stindel, Séverine Ansart. Modeling Early Stages of Bone and Joint Infections Dynamics in Humans: A Multi-Agent, Multi-System Based Model. Frontiers in Molecular Biosciences, 2020, 7, pp.26. 10.3389/fmolb.2020.00026. hal-05247866

HAL Id: hal-05247866

https://hal.univ-brest.fr/hal-05247866v1

Submitted on 10 Sep 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Check for
updates

# Modeling Early Stages of Bone and Joint Infections Dynamics in Humans: A Multi-Agent, Multi-System Based Model

Salma Alsassa¹,²*, Thomas Lefèvre³,⁴, Vincent Laugier², Eric Stindel¹,⁵ and Séverine Ansart¹,⁵

¹ Laboratory of Medical Information Processing (LaTIM - UMR 1101 INSERM), IBRS, Université de Bretagne Occidentale, Department of Medicine, Brest, France, ² Tekliko SARL, Paris, France, ³ Iris UMR 8156 CNRS - U997 Inserm - EHESS - UP 13, Paris, France, ⁴ AP-HP, Jean Verdier Teaching Hospital, Department of Legal and Social Medicine, Bondy, France, ⁵ La Cavale Blanche University Hospital, Infection Diseases Unit, Brest, France

Diagnosis and management of bone and joint infections (BJI) is a challenging task. The high intra and inter patient's variability in terms of clinical presentation makes it impossible to rely on a systematic description or classical statistical analysis for its diagnosis. Advances can be achieved through a better understanding of the system behavior that results from the interactions between the components at a micro-scale level, which is difficult to mastered using traditional methods. Multiple studies from the literature report factors and interactions that affect the dynamics of the BJI system. The objectives of this study were (i) to perform a systematic review to identify relevant interactions between agents (cells, pathogens) and parameters values that characterize agents and interactions, and (ii) to develop a two dimensional computational model of the BJI system based on the results of the systematic review. The model would simulate the behavior resulting from the interactions on the cellular and molecular levels to explore the BJI dynamics, using an agent-based modeling approach. The BJI system's response to different microbial inoculum levels was simulated. The model succeeded in mimicking the dynamics of bacteria, the innate immune cells, and the bone mass during the first stage of infection and for different inoculum levels in a consistent manner. The simulation displayed the destruction in bone tissue as a result of the alteration in bone remodeling process during the infection. The model was used to generate different patterns of system behaviors that could be analyzed in further steps. Simulations results suggested evidence for the existence of latent infections. Finally, we presented a way to analyze and synthesize massive simulated data in a concise and comprehensive manner based on the semi-supervised identification of ordinary differential equations (ODE) systems. It allows to use the known framework for temporal and structural ODE analyses and therefore summarize the whole simulated system dynamical behavior. This first model is intended to be validated by *in vivo* or *in vitro* data and expected to generate hypotheses to be challenged by real data. Step by step, it can be modified and complexified based on the test/validation iteration cycles.

**Keywords:** bone and joint infections, agent-based model, NetLogo, multi-scale, bone remodeling, diagnosis, Staphylococcus aureus

**OPEN ACCESS**

**Edited by:**
Lisa Sedger,
University of Technology Sydney,
Australia

**Reviewed by:**
Ethymia Giannitsioti,
University General Hospital Attikon,
Greece
Haider Abdul-Lateef Mousa,
University of Basrah, Iraq
Bernard Joseph Hudson,
New South Wales Health Pathology,
Australia

*Correspondence:
Salma Alsassa
engsalma80@gmail.com

**Specialty section:**
This article was submitted to Biological Modeling and Simulation, a section of the journal Frontiers in Molecular Biosciences

Received: 10 November 2018

Accepted: 07 February 2020

Published: 12 March 2020

**Citation:**
Alsassa S, Lefèvre T, Laugier V, Stindel E and Ansart S (2020) Modeling Early Stages of Bone and Joint Infections Dynamics in Humans: A Multi-Agent, Multi-System Based Model. *Front. Mol. Biosci.* 7:26. doi: 10.3389/fmolb.2020.00026
# INTRODUCTION

Bone and joint infections (BJI) are challenging pathologies that seriously affect public health, by means of irreversible bone loss, functional impairment and reduced quality of life (Stengel et al., 2001). BJI can result from hematogenous spread of infection and develop within 2 weeks (acute hematogenous BJI) (Carek et al., 2001; Jorge et al., 2010). In adults, BJI develop more frequently from direct inoculation of bacteria secondary to trauma, internal fixation of a fracture, or prostheses placement and progress slowly within months (Smith et al., 2006; Birt et al., 2017). The incidence of BJI is increasing annually in the developed countries due to several reasons such as the rise in arthroplasties, and diabetes (Kurtz et al., 2011; Walter et al., 2012). In 2013, 70 cases per 100,000 of the population were reported as BJI related in France, compared to 54 cases in 2008, along with higher readmission frequency and longer hospital stays. Furthermore, the BJI inpatient mortality rate was reported as 4.6% (Grammatico-Guillon et al., 2012; Laurent et al., 2016). In the USA, Olmsted County, the BJI incidence increased from 11.4 cases per 100,000 of population per year in the period from 1969 to 1979 to 24.4 per year in the period from 2000 to 2009 (Kremers et al., 2015). BJI are also a significant cause of higher healthcare costs each year. The higher cost is due to prolonged and repeated hospitalizations supplemented by complex and long-term treatments, in addition to the need for surgical intervention (Gentry, 1997; Lipsky et al., 2007).

BJI are characterized by a complex biological system on both structural and functional levels. They are defined by multiple variables that have impacts on their pathogenesis (Ciampolini and Harding, 2000). They have a dynamic environment that comprised of different bone matrix qualities, cell types, in addition to the potential presence of biofilms or prostheses (Buckwalter et al., 1996). They encompass complex interactions at both cellular and molecular levels that are difficult to observe and control for using the traditional methods (Lebeaux et al., 2013). Thus, BJI should be better understood through investigating the system behavior that evolves from the interactions between the individual agents. Improving the understanding of the infection dynamics is key to selecting the optimum treatment method.

Numerous animal models have been proposed to study the BJI caused by *Staphylococcus aureus*, the most common cause of this disease (Crémieux and Carbon, 1997; Xing et al., 2012; Johansen et al., 2013). However, these models were limited by the duration of the experiments and the incapacity to observe and control for the interactions between the agents within the process (Mader, 1985; Lebeaux et al., 2013). In addition, the high inter and intra patient's variability of the disease presentation imposes another limitation in front of the classical statistical analysis tools. Thus, the lack of the available evidence and the challenges that accompany the proposed traditional methods necessitate exploring the strength of computational tools to (1) simulate the BJI systems, (2) analyze their progression, and (3) provide insights in selecting evidence-based treatment strategies.

BJI development is the result of the continued interplay between bacteria and two systems within the host: the bone

tissue and the immune defense system. Bacterial invasion of the bone leads to a cascade of adverse changes in bone tissue system components and functions. Specifically, it causes an imbalance in the bone remodeling process, through interactions between the bacteria and bone cells, leading to the activation of bone destruction pathway (Claro et al., 2011; Josse et al., 2015). Concurrently, the presence of bacteria triggers the innate immune cells to signal a response against the invader. The infected bone cells also take part in releasing cytokines that contribute to either directing immune cell reactions or increasing severe inflammatory damage (Bost et al., 1999; Marriott et al., 2005).

Having an integral insight into the disease needs to integrate the various components that manage and direct the dynamic interactions within the BJI process in a computational modeling framework. To date, scarce attempts in the literature investigated the BJI using computational modeling techniques. Pietro and associates (Lio et al., 2012), proposed a computational model to compare the effects of osteomyelitis and osteoporosis on bone remodeling process using the differential equations and probabilistic verification methods based on population-based approach. However, the model emerges from the equation-based approach which lacks the spatial distribution and the micro-interactions of the components, and it ignores the interplay with the immune system.

In this context, the agent-based modeling approach offers a promising framework to translate the system biology and generate a plausible representation of the individual component behavior. ABM approach is a rule-based, bottom-up approach where the system is represented through its entities (agents) that act autonomously in the environment, and the overall behavior of the system results from the accumulation of the local attitudes and interactions of the concerned agents. This approach is based on the principle that mainly local interactions are described. These interactions can take place between agents or between agents and their environment. Global behavior is therefore observed as a phenomenon emerging from the multiple interactions defined in the model. The simulation itself then makes it possible to observe the evolution of the system, by letting it "live" according to the rules and interactions that have been defined on a local scale (An et al., 2009; Shi et al., 2014).

The different features of the ABM approach such as the built-in randomness, spatial architecture, emergent behavior, and taking parallel mechanisms outcomes capability make ABM well suited to represent the biological system and display its hidden and multi-scale behavior (Folcik et al., 2007; Politopoulos, 2007). The ABM models could reproduce the complex behavior of the biological system with heterogeneous components and diverse essential rules, even if they are simple and not complete (Gelfand, 2013; Hammond, 2015). The ABM also has an intuitive paradigm that considers the individual character, decision and propensity. It has the ability to represent the non-linear relationships at the micro-level, which contributes to enhancing the understanding and controlling the different patterns of the system behavior (Kitano, 2002a,b; Ghosh et al., 2011; Hernán, 2015). Using ABM facilitates system simulations to analyze processes that were not
measured in the experimental set-up otherwise (Holcombe et al., 2012; Horn et al., 2012).

In this study, we identified experimental evidence on various aspects of BJI in adults following a literature review and integrated it in a multi-agent, multi-system model based synthesis. We focused on the early stage of the infection to investigate the critical role of innate immune cells in eliminating the bacteria, as well as the loss in bone tissue in this stage. We are focusing on how the cellular mechanisms and interactions influence the BJI development and give rise to the system-level behavior. We developed a two-dimensional agent-based model of BJI that introduces a plausible representation of the system and a simulated experimental environment to reproduce the infection and examine its dynamics as a result of cellular interactions. The BJI system will be modeled through the interaction between the components of three subsystems: the bacteria, the infected bone, and the immune response. Concerning the investigated bacteria, we studied here in this first model S. aureus in general, which is the major bacteria in BJI, notably because it expresses different strategies of persistence such as producing biofilm or hiding intracellularly, and for that it is mainly involved in relapsing infections. The model simulations will be used to test several hypotheses that exist in the literature and to investigate the effect of different parameters changes on the cells dynamics, such as the initial concentration of bacteria, which is a major element in the clinical expression of this infection, and the quality of the immune system during the first stage of infection.

# MATERIALS AND METHODS

In this study, a three-step approach was followed: (i) to conduct a systematic literature search identifying behavior and interactions between agents in both cellular and molecular levels, as well as the values of the parameters involved, (ii) to develop a two-dimensional model of early stages of BJI using the ABM approach, based on the results obtained in the first step, and (iii) to conduct simulations to investigate the outcomes of the infected system under different conditions.

# Literature Review for Agents Interactions and Parameters Values Identification

The first step toward ABM model development was describing the BJI system and its behaviors and characterizing its components and their roles. A systematic literature search was conducted to retrieve the agents' behaviors, rules, and interactions, and to identify the relevant parameters and their ranges of applicable values. By giving some parameters a range of values, we enriched the ABM model with the capacity of simulating variable behaviors of the agents and testing several hypotheses. We used the information from the basic biology to characterize the cellular and molecular attitudes of different components of BJI biological system during the first stage of infection; namely those cover the innate immune response, the bone remodeling cells and signals, and the bacteria development. These biological characteristics will be simplified and formulated as agents' rules. The different parameters were retrieved by

first looking for human studies, but regarding the lack of several information, we extended our search to include general characteristics of the components, animal models of BJI, and mathematical models of the bone remodeling process.

Search was conducted in PubMed and ScienceDirect databases exploring the relevant journal papers from their inception to date, July 2017. Literature search strategy, i.e., search queries and inclusion/exclusion criteria, are described in detail in **Supplementary Material and Supplementary Table S1**.

# The Agent-Based Model of Early Stages Bone and Joint Infections

## The Environment

NetLogo was chosen for implementing the BJI model for its several features such as supporting investigating the agent behaviors over time and space and emerging the agents' heterogeneity within its framework (Wilensky, 1999). Among ABM toolkits, NetLogo introduces an academic and community-supported tool to simulate the compound phenomena using a high-level multi-agent programming language and a robust, simple modeling platform.

The general 2D space of our model was divided into two rectangles of overall dimension 151 × 101 patches (patches represent the grids in the landscape in NetLogo). The first rectangle (101 × 101 patches) was a bone tissue represented as a two-dimensional layer; we assumed this dimension according to 2 mm² of bone tissue. We also assumed that this surface contained one basic multicellular unit (BMU). The second rectangle was the adjacent surface (50 × 101 patches), where the bacteria was initially present (Figure 1A). The virtual time scale was used in the environment “ticks” with each tick mapping to 1 h of the real-time. The model algorithm started with an initial state with all the parameters set to their initial values, followed by a set of rules and functions that were repeated until the end of the simulation.

## Agent Types, Rules, and Parameters

The agents in the model were chosen based on their importance during the first stage of the infection. The bone and immune cell agents interact with the bacteria agents through main interactions while maintaining their functions (Figure 2). Furthermore, the agents in the model and the mediators with their functions and roles are summarized in Tables 1, 2.

### Bone tissue

The bone tissue was represented through the main cells and signals that involved in the bone remodeling process. The agents of bone tissue consisted of Osteoclasts, Osteoblasts, and Osteocytes cells, surrounded by an extracellular matrix (ECM) (the non-cell agents). The interactions between these cells were mediated by the receptor activator of nuclear factor-kappa B ligand (RANKL)/RANK/osteoprotegerin (OPG) signaling system. The individual agents of the three types of bone cells were non-mobile agents and each agent occupied one patch in the grid. The RANKL and OPG signals were represented as mobile agents released by osteoblast cells, whilst the osteoclasts were
<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | Screenshots of the ABM space at three different time steps for inoculum infection state of $(5 \times 10^2$ CFU/mm$^2$). The left rectangle in each sub-figure represents the bacteria population, and the right rectangle represents 2 mm$^2$ of bone tissue where the bone cells, osteoblasts, osteoclasts, and osteocytes, are randomly allocated respecting their percentage and minimum distance. (A) Shows the initial state of the model at $= 0$ h, where the bacteria are randomly distributed in the adjacent surface, the left rectangle, with low presence of immune cells especially macrophages. (B) Shows the state at time $t = 60$ h, where the bacteria entered the bone tissue and started destroying it. (C) Shows the model' state at $t = 150$ h, where the damage happens to the bone tissue, the black patches within bone tissue reflect this destruction, while at the same time the bacteria count was decreased because of engulfing by immune cells.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**FIGURE 2** | Schematic diagram showing agents used in the model, interactions between them, and governing functions for each of them. The oval shapes represent cell agents; hexagons represent signals in the model. The rectangle boxes represent the main functions and roles of each agent. A solid arrow indicates the flow of agent functions, a dotted arrow characterizes stimulation from source to target (destination), while a double-lined arrow reflects the opposite effect (source leads to reduce the destination object).</figcaption></figure>
<table><caption>**TABLE 1** | List of the agents in the bone and joint infections agent-based model, their rules, and behaviors.</caption><thead><tr><td>Agent type</td><td>Agent parameters</td><td>Agent' rules in bond and joint infections ABM</td><td>References</td></tr></thead><tbody><tr><td>Bacteria</td><td>Inoculum size<br>Reproduction rate</td><td>Increase rapidly and spread spatially to invade the bone tissue, stimulate releasing RANKL and activating OC, stimulate immune defense and are engulfed by them, stimulate OB death</td><td>Fullilove et al., 2000; Claro et al., 2011; Olson and Horswill, 2013</td></tr><tr><td>Neutrophils (PMN)</td><td>Count, Lifespan<br>Reproduction rate</td><td>Undergo reproduction and death function, recruited due to the presence of bacteria and try to ingest them, recruit MDM</td><td>Cassatella, 1995; Kaplanski et al., 2003; Bekkering, 2013; Boff et al., 2018</td></tr><tr><td>Macrophages (MA)</td><td>Count, Lifespan<br>Reproduction rate</td><td>Undertake reproduction and death function, stimulated by the presence of bacteria and attack them, regulate macrophages and MDM recruitment through TGF-beta, stimulating neutrophils</td><td>Kim and Luster, 2015; Kaufmann and Dorhoi, 2016; Prame Kumar et al., 2018</td></tr><tr><td>Monocytes (MDM)</td><td>Count, Lifespan<br>Reproduction rate</td><td>Undergo reproduction and death function, stimulated by the presence of bacteria, PMN, and MA after T hours, phagocytosis the bacteria, release TGF-beta to regulate macrophages and MDM recruitment</td><td>Deshmane et al., 2009; Lee et al., 2011; Gong et al., 2012; Ginhoux and Jung, 2014; Patel et al., 2017</td></tr><tr><td>Osteoblasts (OB)</td><td>Count, Lifespan<br>Reproduction rate</td><td>Go through reproduction and death cycle, spatial localization, releasing RANKL and OPG, take a role in bone remodeling process: form new osteocytes</td><td>Jilka et al., 1998; Boyce and Xing, 2008; Shapiro, 2008</td></tr><tr><td>Osteoclasts (OC)</td><td>Count, Lifespan<br>Reproduction rate</td><td>Go through reproduction and death cycle, spatial localization, bind with RANKL to be activated, take a role in bone remodeling process: destroying osteocytes</td><td>Parfitt, 1994; Blair, 1998; Eriksen, 2010</td></tr><tr><td>Osteocytes (OS)</td><td>Count, percentage</td><td>Form bone osteocytes cells network with respecting the minimum distance between them, derived from mature osteoblasts, destroyed by active osteoclasts</td><td>Franz-Odendaal et al., 2006; Bonewald, 2007; Schneider et al., 2010</td></tr></tbody></table>

PMN, Neutrophils; MA, Macrophages; MDM, Monocyte-derived Macrophages; TGF-β, Transforming growth factor-beta; RANKL, Receptor activator of nuclear factor-kappa-B ligand; OPG, Osteoprotegerin.

<table><caption>**TABLE 2** | List of the mediators and their effects that are represented in the bone and joint infections ABM.</caption><thead><tr><td>Mediator variable</td><td>Mediator parameter</td><td>Source</td><td>Role in BJL agent-based model</td><td>References</td></tr></thead><tbody><tr><td>RANKL</td><td>Concentration</td><td>Osteoblasts</td><td>Diffusion, activate osteoclasts by binding to them or inhibit their activation by binding with OPG</td><td>Anandarajah, 2009</td></tr><tr><td>OPG</td><td>Concentration</td><td>Osteoblasts</td><td>Diffusion, bind with RANKL to inhibit activating osteoclasts</td><td>Boyce and Xing, 2008</td></tr><tr><td>TGF-beta</td><td>Concentration</td><td>Macrophages<br>Monocytes</td><td>Released by both monocyte and macrophages to increase monocytes recruitment and decrease macrophage recruitment</td><td>Ashcroft, 1999</td></tr><tr><td>MCP-1</td><td>Concentration</td><td>Neutrophils</td><td>Released by neutrophils to stimulate monocytes recruitment</td><td>Gopalakrishnan et al., 2013; Tecchio et al., 2014</td></tr><tr><td>TNF</td><td>Concentration</td><td>Macrophages</td><td>Released by macrophages to enlist neutrophils to the site</td><td>Prame Kumar et al., 2018</td></tr></tbody></table>

BJL, bone and joint infections; TGF-β, Transforming growth factor-beta; RANKL, Receptor activator of nuclear factor-kappa-B ligand; OPG, Osteoprotegerin; MCP-1, Monocyte chemoattractant protein-1; TNF, Tumor necrosis factor.

modeled to react as RANK receptors. The osteoclast agents were modeled to perform bone resorption through destroying existing osteocyte agents, after being activated by binding to RANKL. On the other hand, the osteoblast agents were modeled to move toward building new osteocyte agents within the cell-matrix and to control the osteoclast agents' activation through RANKL/OPG signals. While RANKL was modeled to increase osteoclasts activation, OPG, in turn, is modeled to inhibit this activation by binding to RANKL. In a healthy state, the resorption and formation of new osteocytes are balanced. However, this process is modeled to be altered during the infection by increasing the RANKL concentration and subsequently the osteoclasts activity toward increasing bone destruction. Bone cells were modeled to have the following behavior rules and functions.

### Osteoclasts.

1. During initialization, they are distributed depending on their count in one BMU (Komarova et al., 2003; Paoletti et al., 2012).

2. They go through apoptosis when their ages, which are increased by each time step, reach the lifespan parameter (Manolagas, 2000; Ryser et al., 2009; Lio et al., 2012).

3. New osteoclasts are reproduced every day depending on their reproduction rate parameter.

4. The osteoclasts have RANK receptors on their surfaces (Roodman, 1999), which cause osteoclasts activation through binding to RANKL molecules.

5. Activated osteoclasts move toward destroying neighboring osteocytes cells (Figure 2).

### Osteoblasts.

1. Their initial distribution is based on their count in one BMU (Komarova et al., 2003; Paoletti et al., 2012; Florencio-Silva et al., 2015).

2. When their ages reach their lifespan parameter, they change their status to new osteocytes (Dallas and Bonewald, 2010).
3. New osteoblasts are reproduced every day depending on their reproduction rate parameter.

4. Osteoblasts agents are considered as the source of releasing RANKL. Releasing RANKL is depending on both the reproduction rate parameter and the existence of the bacteria. RANKL modeled to bind to RANK and activate osteoclasts (Figure 2). The releasing process was verified using diffuse instruction in NetLogo.

5. OPG molecules are also released by osteoblasts agents using their reproduction rate and modeled to bind to RANKL in order to inhibit osteoclasts activation.

6. Osteoblasts agents migrate toward forming new osteocytes, after searching for a place where no osteocyte occupies a nearby location to maintain a minimum distance between two osteocytes (Repp et al., 2017).

### Osteocytes.

1. At the initial state, they are distributed in the bone tissue space taking into account a minimum distance between them forming a network of osteocytes cells within the extracellular matrix of bone.

2. Because of their long lifespan that could extend for decades (Franz-Odendaal et al., 2006; Dallas et al., 2013), they were not modeled to undergo reproduction and death phase.

3. Instead, new osteocytes originated from mature osteoblasts and were destroyed by activated osteoclasts.

### Immune system

It was modeled with three types of innate immune cells namely macrophages, neutrophils, and monocyte-derived macrophages (MDM), which represent the first actors against the bacteria and they respond to the infection in sequential roles (Kaplanski et al., 2003; Strydom and Rankin, 2013; Kim and Luster, 2015). These immune cells were represented as mobile agents “turtles” in the model. The response sequence is modeled on two steps. First, the macrophage and neutrophil agents are modeled to react directly against the bacteria and phagocytose them. In addition, their stimulation is increased by the proliferation of bacteria and regulated by the tumor necrosis factor TNF and transforming growth factor-beta TGF-β cytokines. The second step of defense is modeled to take place if both of macrophages and neutrophils cells could not succeed in eliminating the bacteria by the first 48 h, through recruiting MDM cells to the site of infection. The MDM cells were modeled to be stimulated by the Monocyte chemoattractant protein-1 (MCP-1) cytokine. MDM population is also regulated by the TGF-β signal through a positive feedback loop. TGF-β was represented as mobile agents in the model, while MCP-1 and TNF were represented as variables associated with the infected patches. The functions and rules behind the innate immune cells' work follow.

### Macrophages.

1. At initialization, macrophages agents are distributed randomly in the model space depending on their initial number, and they are stimulated by the presence of the bacteria.

2. Macrophages go through apoptosis depending on their lifespan parameter, while their ages increase at each time step.

3. New macrophages are reproduced depending on their reproduction rate. The reproduction process follows the uniform distribution function in which the probability of reproducing new macrophages in one time interval is the same for the whole reproduction time.

4. The macrophages agents move around the space and are attracted to the presence of bacteria to phagocytosis them. When one macrophage engulfs a bacterial agent, it dies.

5. The macrophages reproduce TGF-β cytokines that have a double role. They decrease the reproduction of macrophages through a negative feedback loop, while on the second hand, they increase the monocyte production.

6. The macrophages also increase the inflammation by releasing the TNF cytokine as an associated variable that used to interact with neutrophils agents to increase their recruitments.

### Neutrophils.

1. Since they are recruited to the infection site within 4–8 h (Strydom and Rankin, 2013), they were modeled to react against the invader with no delay. Neutrophils agents are created with the initialization of the model and randomly distributed in the modeled area with an initial concentration.

2. Neutrophils go through apoptosis depending on their lifespan parameter, while their ages increase at each time step.

3. New neutrophils agents are produced each day depending on their reproduction rate. Agents' reproduction is also equally distributed over the time of reproduction.

4. Neutrophils agents follow a random walk toward phagocytosis the bacteria. They die after engulfing the bacteria.

5. They are stimulated by the presence of bacteria and by the pro-inflammatory cytokine TNF that are released by the macrophages. In addition, they are considered as the source of MCP-1, through which the monocytes are recruited.

### Monocyte-derived-macrophages (MDM).

1. They are recruited to the site 48 h post the bacterial invasion, stimulated by macrophages and neutrophils through MCP-1 and TGF-β cytokines.

2. They have a life cycle that goes through increasing the age and ends by the apoptosis depending on the lifespan parameter (Whitelaw and Bell, 1966; Ginhoux and Jung, 2014; Italiani and Boraschi, 2014).

3. Each day, new MDM agents are modeled to be produced following the uniform distribution function.

4. They are modeled to move toward phagocytosis the bacteria, and they die after engulfing one.

5. They were also considered as a source of TGF-beta signal, through which their count was increased through a positive feedback loop.
The reproduction rate of each of the immune cells was estimated using the method used in the mathematical model of Smith et al. (2011). We assumed that the steady-state $N$ is given by $N = s/d$, where $s$ is the reproduction rate of cells per day and $d$ is the clearance rate per day. The clearance rate of macrophages, which is taken from their identified lifespan, is in the range $d = 7-14 \times 10^{-2} \text{ day}^{-1}$, implying $s = 28-115 \text{ cell/day}$, for a steady-state 200-800 cell/mm³ (Wallace et al., 1993). In the same way, the reproduction rate for each of neutrophils and monocytes was estimated for steady values taken from an animal model of bone infections (Corrado et al., 2016). The reproduction rate of neutrophils was calculated as 120-700 cell/day for steady-state 250-700 cell/mm³ and $d = 0.2-1 \text{ day}^{-1}$. When the monocytes are activated, they will reproduce at a rate of 4-70 cell/day that estimated for $d = 0.2-1 \text{ day}^{-1}$ and steady-state 20-70 cell/mm³.

## Bacteria

The bacteria in the model were modeled to be affected by several factors; some of them were related to the bacteria themselves, such as their initial inoculum size and reproduction rate. The other factors were those related to the immune cells' defense and ability to eliminate the bacteria. In their turn, the bacteria dynamics affected the bone tissue health represented by the number of osteocytes or ECM agents. The bacteria were represented in the model as turtles with each individual mapped to one patch, and they propagate depending on their reproduction rate parameter (Table 3). The reproduction rate of bacteria differs in the controlled culture from it in the human body. This rate is difficult to be determined in the human body because of the different factors that have impacts on the bacterial growth in the infected bone such as the location, vascularization, pH, nutrition, and type of prosthesis, if exist. The bacteria characteristics, such as SCV, also have an impact on the growth rate toward decreasing it (Bui et al., 2015). The reproduction rate in the model was identified by the range of 1-24 h to cover several proposed values (Fux et al., 2005; Anwar et al., 2007). The modeled bacterial behaviors rules follow:

1. Bacteria agents are randomly distributed depending on their initial inoculum size.

2. Their life cycle goes through the reproduction phase depending on their reproduction rate parameter.

3. The bacteria spread and move toward invading the bone tissue.

4. The bacteria go through death rate which symbolizes the run out of their resources of survival.

5. The existence of bacteria stimulates the immune cells to react against the invasion starting by macrophages and neutrophils agents.

## Process Overview and Events Sequence

The simulation was launched using a function that mainly sets all the variables to their initial values, in addition to carrying out the following tasks. First, it established a basic grid of patches (151 × 101 patches). Second, it initiated the bone cells patches according to their percentage of bone cells count and allocated them randomly. Third, it created

numbers of macrophage and neutrophil agents corresponding to their initial value and distributed them randomly. Fourth, it created a number of bacteria agents (turtles) respecting their chosen initial value, and randomly allocated each cell to one patch on the adjacent part of the grid. Finally, it created numbers of signaling agents (TGF-β, RANKL, OPG) corresponding to their initial concentration level, and distributed them randomly.

After initialization, the simulation started looking through the set of functions and rules over time and were repeated up to the end of simulation either by eliminating the bacteria or reaching the determined simulation time (300 h). After the initialization of bacteria, proliferation was modeled using binary fission depending on varied reproduction time. The fission function was modeled to follow the normal distribution over the time of reproduction (Table 3). For one bacteria to undergo fission, it must search for one empty neighboring patch of the bone tissue; otherwise, no fission would take place. At each time step, the agents updated their location and moved upon their functions and rules, the bacteria moved toward invading the bone tissue, the immune cells moved toward engulfing the bacteria, and the osteoclasts and osteoblasts updated their location toward resorbing and forming new osteocytes.

The population and distribution of each type of agent were calculated at each time step, if some agents reached their lifespan, they went under apoptosis, otherwise, they increased their age. At the time of reproduction, new agents were created upon their reproduction rate and allocation conditions. During their life, each type of agents performed its rules and functions that described in section Agent Types, Rules, and Parameters.

During the simulation and at each time step, the populations of a pre-defined set of agents were saved to ".csv" file. These output files were identified by the initial conditions of the simulation. Multiple runs of the same simulation were saved to the same output file in order to analyze the model outcomes for the same initial conditions. The progress of the infection over time and space during the simulation is shown in the simulator interface (Figure 1). Changes in the agents' counts were monitored through time graphs for each type of agent in the NetLogo interface.

## Simulated Experiments

The goals of the simulation were to illustrate the relationships between variation in the bone cell population and the evolution of infectious processes. The goal was also to investigate the efficiency of the innate immune system in defending the bacteria during the first stage of the BJI for several scenarios without treatment intervention. Further, the simulations were conducted to detect sensitive parameters that affected the dynamics of the system.

To verify the model stability and performance, and to investigate the effect of bacterial inoculum size on the system dynamics, the initial number of bacteria were adjusted to three different values (5, 5 × 10, 5 × 10² CFU/mm²). These values were assumed to represent three different infected inoculum states – low (5 CFU/mm²), medium (5 × 10 CFU/mm²) and
<table><caption>TABLE 3 | The parameters in the model and their values or ranges used.</caption><thead><tr><th>Parameter</th><th>Range in literature</th><th>References</th><th>Type of study</th><th>Range in the model</th><th>Step size</th><th>Simulation value</th></tr></thead><tbody><tr><td>Bacteria production-rate</td><td>1–24 h</td><td>Fux et al., 2005; Anwar et al., 2007</td><td>In vitro</td><td>[1–24] hour</td><td>1 h</td><td>12 h</td></tr><tr><td>Bacteria inoculum size</td><td>0–500 CFU/mm³</td><td>Assumed</td><td></td><td>[0–500] CFU/mm²</td><td>10 CFU/mm²</td><td>5, 50, 500 CFU/mm²</td></tr><tr><td>Osteocytes initial number</td><td>500–900 cell/mm²</td><td>Vashishth et al., 2000; Goggin et al., 2016</td><td>In vivo, Human study</td><td>1500–2000 cells</td><td>–</td><td>1790 cells</td></tr><tr><td>Osteoblasts production-rate</td><td>4 cell/day</td><td>Parfitt, 1994; Komarova et al., 2003</td><td>In vivo, human study</td><td>[1–10] cell/day</td><td>1 cell/day</td><td>4 cell/day</td></tr><tr><td>Osteoblasts lifespan</td><td>3 months</td><td>Manolagas, 2000</td><td>Human study</td><td>[10–90] day</td><td>5 days</td><td>50 days</td></tr><tr><td>Osteoblasts initial number</td><td>800–2000 cells/BMU</td><td>Parfitt, 1994; Komarova et al., 2003</td><td>In vivo, human study</td><td>800–2000 cells</td><td>–</td><td>1000 cells</td></tr><tr><td>Osteoclasts production-rate</td><td>3 cell/day</td><td>Parfitt, 1994; Komarova et al., 2003</td><td>In vivo, human study</td><td>[1–5] cell/day</td><td>1 cell/day</td><td>3 cell/day</td></tr><tr><td>Osteoclasts lifespan</td><td>2 weeks</td><td>Manolagas, 2000</td><td>Human study</td><td>[1–14] day</td><td>1 day</td><td>7 days</td></tr><tr><td>Osteoclasts initial number</td><td>5–20 cells/BMU</td><td>Parfitt, 1994; Komarova et al., 2003</td><td>In vivo, human study</td><td>5–20 cells</td><td>–</td><td>8 cells</td></tr><tr><td>RANKL concentration</td><td>10<sup>-6</sup> mol/cell/day</td><td>Ryser et al., 2009</td><td>Estimated according to in vivo observation</td><td>1 µmol/cell/day</td><td>1 µmol/cell/day</td><td>1 µmol/cell/day</td></tr><tr><td>OPG concentration</td><td>3.10<sup>-6</sup> mol/cell/day</td><td>Ryser et al., 2009</td><td>Estimated according to in vivo observation</td><td>3 µmol/cell/day</td><td>1 µmol/cell/day</td><td>3 µmol/cell/day</td></tr><tr><td>TGF-β concentration</td><td>150–500 pg/ml</td><td>Knapp et al., 1998</td><td>In vivo, human study</td><td>1 × 10<sup>-3</sup> pg/cell/day</td><td>–</td><td>1 × 10<sup>-3</sup> pg/cell/day</td></tr><tr><td>TNF concentration</td><td>0–1000 pg/ml</td><td>Corrado et al., 2016</td><td>In vivo, animal study (murine)</td><td>1 × 10<sup>-3</sup> pg/cell/day</td><td>–</td><td>1 × 10<sup>-3</sup> pg/cell/day</td></tr><tr><td>MCP-1 concentration</td><td>0–2000 pg/ml</td><td>Corrado et al., 2016</td><td>In vivo, animal study (murine)</td><td>1 × 10<sup>-3</sup> pg/cell/day</td><td>–</td><td>1 × 10<sup>-3</sup> pg/cell/day</td></tr><tr><td>Neutrophil reproduction-rate</td><td>–</td><td>Estimated</td><td></td><td>[120–700] cell/hour</td><td>50 cells</td><td>550 cell/day</td></tr><tr><td>Neutrophil lifespan tissue</td><td>24–120 h</td><td>Rankin, 2010; Bekkering, 2013</td><td>In vivo, human study, mice study</td><td>[24–120] hour</td><td>6 h</td><td>60 h</td></tr><tr><td>Monocyte lifespan</td><td>24–120 h</td><td>Whitelaw and Bell, 1966; Ginhoux and Jung, 2014; Italiani and Boraschi, 2014; Patel et al., 2017</td><td>In vivo, human study, mice study</td><td>[24–120] hour</td><td>5 h</td><td>60 h</td></tr><tr><td>Monocyte reproduction rate</td><td>–</td><td>Estimated</td><td></td><td>[4–70] cell/day</td><td>50 cell/day</td><td>150 cell/day</td></tr><tr><td>Macrophage lifespan</td><td>1–14 days</td><td>Parwaresch and Wacker, 1984</td><td>In vivo, rat study</td><td>[24–300] hours</td><td>6 h</td><td>24 h</td></tr><tr><td>Macrophage reproduction-rate</td><td>–</td><td>Estimated</td><td></td><td>[28–115] cell/day</td><td>10 cell/day</td><td>550 cell/day</td></tr></tbody></table>

BMU, Basic Multicellular Unit.

high (5 × 10² CFU/mm²). We ran the simulations for n = 100 iterations for each of the infected inoculum state, while all other changeable parameters were set to the median of their range (Table 3). For each iteration, the dynamics of the population for bacteria, osteocytes, and neutrophils agents for a time duration of 300 h (t = 300 ticks) was tracked. The generated data at each time click (1 h) were saved to an output data (.csv) file for analysis. Subsequently, the mean and standard deviation for the population at each time step for each type of these agents was quantified. Further, the relation between two agents over time under the same initial conditions was analyzed. We benefitted from the 3D surface graphs and the generated data from the previous step to analyze the relationships between each of bacteria and neutrophils populations over time, bacteria and

osteocytes populations over time, and neutrophils and osteocytes populations over time.

## RESULTS

### Characterizing the Agent Interactions and Parameters Values

The systematic search process identified 42 articles for bone cell characterizations, 48 articles for immune cell characterization, 12 articles for bacteria characterization, and 29 articles for interactions (Supplementary Figures S1–S4). Following the systematic literature search, the agent interactions and the parameter ranges were identified from the retrieved articles.
Flowcharts of the review steps along with the list of all articles used for the study are reported in **Supplementary Figures S1-S4** and **Supplementary Table S2**.

## The Developed Agent-Based Model of BJI

A representative example of the model space was taken during one simulation at three different times (Figure 1). The initial state of the model was represented at the beginning of the simulation (0 h), for initial inoculum infection state of $(5 \times 10^2 \text{ CFU/mm}^2)$, where this initial number of bacteria were randomly distributed in the adjacent surface. It also demonstrated the low presence of immune cells macrophages and neutrophils, and the initial state of bone tissue before being damaged (Figure 1A). The progress of the infection during the simulation (60 h) later showed that the bacteria propagated toward entering the bone tissue and destroying it. The monocytes also had been activated in this step besides the other immune cells (Figure 1B). The simulation state at (150 h) showed the system after 6 days, where the damage happened to the bone tissue, and the bacteria count decreased because of engulfing by the immune cells (Figure 1C). This model is 2D what makes it suitable to represent trabeculae bone, but it could be also considered as a cross-section of cortical bone.

## Simulated Experiments Output

The mean and the standard deviation of the dynamics of the population for 100 iterations under different infection inoculum states (low, medium, high) revealed interesting behavior for each of the bacterial, neutrophil (PMN), and osteocyte population (Figure 3). The bacterial population inclined toward the same steady non-null counts, regardless of the inoculum state (Figures 3A1–A3). In addition, the population intensity of the first few days was proportional to inoculum size. It was also observed that the behavior of the bacteria varied compared to the mean behavior with small variance magnitude even if the output illustrated (by the effect of time characteristic) a high-frequency oscillation in the bacteria population. At the same time, it showed low fluctuation in the context of the general trend of the bacterial population (Figures 3A1–A3). It was also observed that the bacterial population faced quasi-extinction, followed by regrowth for the inoculum in the medium and high infection state for all 100 iterations.

On the other hand, the PMN dynamics followed the evolution of the bacteria with a slight delay, reaching a non-zero stable level on the 12th day. In addition, it illustrated that the PMN population was not subject to rapid variations represented by the low-frequency oscillations, which are comparatively smoother than those introduced by the bacteria (Figures 3B1–B3). Asymptotic behavior of PMN population for the three inoculum states proposed a non-null mean value with fixed frequency oscillations and decreasing magnitude with time (Figures 3B1–B3).

Concerning the bone tissue loss in each state, the osteocyte population dynamics was similar in mean population intensity, while it differed from the mean behavior with important variance magnitudes whatever the inoculum. The bone cell population

predicted a trend of unexpected high-frequency oscillations after $t = 200 \text{ h}$ (Figures 3C1–C3).

Comparisons of system response for different inoculum states and over time for the same type of agent population suggested that the high inoculum infection state was associated with close to elimination of bacteria population and full elimination of PMN population. It was also associated with the highest population counts for both agent types in the first few days of the infection (Figure 3). The results presented that bacteria and PMN tended toward a steady asymptotic non-null state, regardless of the inoculum infection state. With regards to bone tissue damage, the first response stage had a similar degradation phase of the bone infected site represented by decreasing count of osteocytes for each inoculum state, where it was noted that osteocytes passed by a common minimum level by the 7th-8th day for all three inoculum states, representing 2% of the loss in the infected site mass (Figure 3). During the recovery phase, it was observed that the smallest inoculum infection state caused the most severe and relatively stable loss on bone cells by 12th day, the medium inoculum infection state showed better progressive recovery, and the high inoculum infection state displayed an intermediate recovery, with sub-optimum recovery at 300 h. The displayed loss of bone tissue compared to baseline for each inoculum infection state was 1.4, 1.2, and 1%, respectively (Figure 3).

The 3D representation of the relation between two types of agents over time suggested minimum levels of bacteria and PMN population count around the fourth day of infection, while it suggested a delay in minimum levels of population in osteocytes count with regard to the minimum levels of bacteria and PMN (Figure 4).

## DISCUSSION

The incidence rate of post-operative/post-traumatic infections is showing an important annual increase due to the rising in arthroplasties procedures associated with mounting risk factors such as diabetes and peripheral vascular disease, in addition to the variation in the population's age structure (Walter et al., 2012). BJI are complex and characterized by different patterns of bacteria progression and different hypotheses in the literature. The lack of strong evidence is exacerbated by indefinite possible dynamic processes and interactions at cellular and molecular levels. This infection is characterized by high divergence during inter and intra infected cases, and by multi pre-defined variables, making it difficult to rely on classical statistical models. By integrating the ABM modeling technique along with the knowledge from the literature review, this study introduced an innovative *in silico* experimental environment to explain the hypotheses and knowledge in the literature and investigate the role and dynamics of other interactions that could be proposed. The developed model offers a means to test the impact of several factors on the infected tissue. Through this model, this study illustrated the role of bacterial virulence and host defense state that identify the consequences of bacterial invasion to the bone. Although the simulations are at the preliminary levels with regards to their accuracy and validity, several interesting
<figure><img src="image_5.png" /><figcaption>**FIGURE 3** | The mean and standard deviation (SD) for 100 iterations for agent populations over time, *t* = 300 h, at three inoculum infection states of bacteria (5, 50, 500 CFU/mm²). (A1–A3) Show the mean and SD of the bacteria population over time. (B1–B3) Show the mean and SD corresponding to neutrophils population. (C1–C3) Represent mean and SD for osteocytes population. First, second and third columns represent the results for inoculum infection states of low (5 CFU/mm²), medium (5 × 10 CFU/mm²), and high (5 × 10² CFU/mm²) bacteria respectively.</figcaption></figure>

outcomes were observed that would certainly provide insights
into the BJI dynamics.

Observation

For the proposed inoculum infection states, on the 12th day,
the observed bacterial population for the inoculum of the
order of 100 was stable and not extinct, which could also
be considered as an indication for latent infection. Since the
PMN populations followed, in the same manner, the bacteria
populations but in a slight delay to reach the stable non-null
levels on the 12th day, it could also be used as an indication
of the latent infection. If so, it should be validated in terms of
PMN count in the biological laboratory test. In addition, the
relationship between the osteocytes and the inoculum infection
state raises the question of whether it states *a restitutio ad integrum* for the medium infection state that introduced the best
progressive recovery.

Further, the observed stable non-extinct level of bacteria
population 10 days after the infection showed the inability of the

innate immune cells to completely eliminate the bacteria. This
points out the imperative role of further defense methods: the
adaptive immune response, which initiated by the 4th–7th day
post the infection, and the need for therapeutic intervention.

In terms of population fluctuations (intra-simulation), the
observed small variance within the bacterial population dynamics
highlighted their stable behavior even with the rapid, strong
oscillation which could be an artifact of the internal clock
management in NetLogo algorithm. The fluctuations observed in
the PMN population had a less important effect compared to the
bacterial state, which reflects a spatial impact. In fact, greater the
probability that a PMN encountered a bacterium, more likely it
would have reacted, which ultimately reduced the fluctuations.

On the other hand, the predicted high-frequency oscillating
trend of the osteocytes population after *t* = 200 h was unexpected.
It was expected to display a rapid decreasing phase, and then
a slower increasing phase. This considerable fluctuations in the
osteocytes population appeared interesting as it could partially
explain the inter-individual variability. Further, these fluctuations
<figure><img src="image_6.png" /><figcaption>**FIGURE 4** | 3D surface graphs to analyze the relationships between two types of agents over time at inoculum infection states of bacteria (50 CFU/mm²). In the top, the graph shows the relationship between bacteria vs. osteocytes (OS) population over time for the three initial inoculum values. The second row corresponds to the osteocytes (OS) vs. neutrophils (PMN) population over time. The third row shows the relationship of neutrophils (PMN) vs. bacteria populations over time.</figcaption></figure>
increased in intensity over time, raising the question if it could
be explained as a pre-chaotic behavior “positive feedback.” It
was also noted that the PMN and the bacteria populations reach
their minimum levels simultaneously before rising again. This
weakness in the innate immune response could provide a better
possible window to start the therapeutic treatment.

Strengths

The BJI system model developed in this study is the first
to incorporate experimental data with the ABM modeling
approach and extract information in the form of a dynamic
system and it lays the foundations for an in-depth and detailed
BJI simulation system. The interactions between the agents
and signals provide a comprehensive ability to analyze the
system considering the spatial characteristics and the inter-
agent variability. Additionally, the observed results corroborate
with clinically and microbiologically available pathophysiological
pathways and they are generally in line with expectations, as that
the investigated populations move toward equilibrium (Wagner
et al., 2003). This model has also shown its predictive ability for
the evolution of bone mass with respect to bacterial inoculum
states and time. The proposed model framework introduces a
flexible and interactive virtual laboratory to test and explain
several existing hypotheses or even to design new ones.

Limitations

One of the major limitations of the ABM approach is that it
generally cannot be released simply, unless the entire model is
given. Evidence is also difficult to obtain, and implementation
details can cause many problems. For that, we had to start with
a less complex model (at the same time integral one) and to build
later upon it. The developed model presents certain limitations,
and all the results should be evaluated within the framework of
these limitations. The first limitation is model validation due to
the lack of experimental data, which might be acquired at a later
stage (e.g., the bone mass). The 2D representation represents a
limitation with regards to the ability of agents to interact in a 3D
space. However, we adopted the simplicity in building this first
model aiming at developing a feasible modeling framework and
understanding the complex integration of available physiological
data with the ABM modeling framework. Future work is aimed
at enhancing the model by a more realistic architecture, tissue
specification, and patient-based data. The current model lacks
the second stage of the immune response, the adaptive immune
response, which is necessary to investigate the further progression
of the infection. The bacteria biofilms, which are lacking in the
current model, also play an important role in identifying the
behavior of the bacteria itself and the system response due to
their resistance to immune defense and antibiotic agents. Since
biofilms are common during this infection, it is worthwhile to
model them in the future. In this model we were limited to
S. aureus, other pathogens will be investigated in the future. In
this first model, the BJI was modeled without prostheses, which
will be added in the future. Other agents including bacterial
survival factors have a significant contribution in introducing
an additional approximation to mimic the real system. The
adaptability of the present model depends on the type of

modifications that are foreseen. All modifications needed would
be investigated through the same process we described in this
study: first, identify if (i) new agents (e.g., cells) should be
integrated and (ii) already modeled interactions are sufficient
or should be modified or completed by new ones; second, list
all relevant parameters and parameters values related to these
modifications through literature search; finally, model agents,
messengers and implement parameters values in the new model.
This process can be applied if we want to introduce variants or
new pathogens that may interact with each other and with the
existing agents or if we want to model cancellous bone rather
than cortical bone or if children and adults BJIs differ from
critical aspects.

Vision

This work has drawn attention toward the incorporation
between modeling power and current literature using the ABM
modeling framework of infection dynamics. It also highlights the
progression of the simulation environment that generates data
that would be used in the extraction and synthesis of the model
in the form of differential equation systems. This model could in
fine be integrated within a global predictive approach to propose
more personalized therapeutic treatments. The future global
approach would combine microbiological and imaging data of
the patient to extract infection parameters (bacteria, immune
system) and patient-specific bone morphology. As further data
are made available, the model would be refined to better
estimate an in-depth comprehension of BJI's pathophysiology.
Furthermore, studying the dynamics of the BJI including critical
parameters such as the agent interactions and cross-talk signals
could make it possible to test new or alternative therapeutic
pathways for their efficacy.

AUTHOR CONTRIBUTIONS

SAI drafted the manuscript and was responsible for applying the
computational algorithm, implementing the ABM model, doing
the simulations, and processing the results. TL, SAn, and ES
were responsible for revising the manuscript, for the design and
interpretation of the simulation, and ABM results. TL, VL, SAn,
and ES supervised the work.

FUNDING

This work was supported by the Tekliko SARL and the
Laboratory of Medical Information Processing (LaTIM -
UMR 1101 INSERM).

SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online
at: https://www.frontiersin.org/articles/10.3389/fmolb.2020.
00026/full#supplementary-material
REFERENCES

An, G., Mi, Q., Dutta-Moscato, J., and Vodovotz, Y. (2009). Agent-based models in translational systems biology. *Wiley Interdiscip. Rev. Syst. Biol. Med.* 1, 159–171. doi: 10.1002/wsbm.45

Anandarajah, A. P. (2009). Role of RANKL in bone diseases. *Trends Endocrinol. Metab.* 20, 88–94. doi: 10.1016/j tem.2008.10.007

Anwar, H. A., Aldam, C. H., Visuvanathan, S., and Hart, A. J. (2007). The effect of metal ions in solution on bacterial growth compared with wear particles from hip replacements. *J. Bone Joint Surg. Br.* 89, 1655–1659. doi: 10.1302/0301-620X.89B12.19714

Ashcroft, G. S. (1999). Bidirectional regulation of macrophage function by TGF-β. *Microbes Infect.* 1, 1275–1282. doi: 10.1016/s1286-4579(99)00257-9

Bekkering, S. (2013). Another look at the life of a neutrophil. *World J. Hematol.* 2, 44–58. doi: 10.5315/wjh.v2.i2.44

Birt, M. C., Anderson, D. W., Bruce Toby, E., and Wang, J. (2017). Osteomyelitis: recent advances in pathophysiology and therapeutic strategies. *J. Orthop.* 14, 45–52. doi: 10.1016/j.jor.2016.10.004

Blair, H. C. (1998). How the osteoclast degrades bone. *Bioessays* 20, 837–846. doi: 10.1002/(sici)1521-1878(199810)20:10<837::aid-bies9>3.0.co;2-d

Boff, D., Crijns, H., Teixeira, M. M., Amaral, F. A., and Proost, P. (2018). Neutrophils: beneficial and harmful cells in septic arthritis. *Int. J. Mol. Sci.* 19:468. doi: 10.3390/ijms19020468

Bonewald, L. F. (2007). Osteocytes as dynamic multifunctional cells. *Ann. N. Y. Acad. Sci.* 1116, 281–290. doi: 10.1196/annals.1402.018

Bost, K. L., Ramp, W. K., Nicholson, N. C., Bento, J. L., Marriott, I., and Hudson, M. C. (1999). Staphylococcus aureus infection of mouse or human osteoblasts induces high levels of interleukin-6 and interleukin-12 production. *J. Infect. Dis.* 180, 1912–1920. doi: 10.1086/315138

Boyce, B. F., and Xing, L. (2008). Functions of RANKL/RANK/OPG in bone modeling and remodeling. *Arch. Biochem. Biophys.* 473, 139–146. doi: 10.1016/j.abb.2008.03.018

Buckwalter, J. A., Glimcher, M. J., Cooper, R. R., and Recker, R. (1996). Bone biology. I: structure, blood supply, cells, matrix, and mineralization. *Instr. Course Lect.* 45, 371–386.

Bui, L. M. G., Hoffmann, P., Turnidge, J. D., Zilm, P. S., and Kidd, S. P. (2015). Prolonged growth of a clinical staphylococcus aureus strain selects for a stable small-colony-variant cell type. *Infect. Immun.* 83, 470–481. doi: 10.1128/IAI.02702-14

Carek, P. J., Dickerson, L. M., and Sackier, J. M. (2001). Diagnosis and management of osteomyelitis. *Am. Fam. Physician* 63, 2413–2420.

Cassatella, M. A. (1995). The production of cytokines by polymorphonuclear neutrophils. *Immunol. Today* 16, 21–26. doi: 10.1016/0167-5699(95)80066-2

Ciampolini, J., and Harding, K. G. (2000). Pathophysiology of chronic bacterial osteomyelitis. Why do antibiotics fail so often? *Postgrad. Med.* J. 76, 479–483. doi: 10.1136/pmj.76.898.479

Claro, T., Widaa, A., O'Seaghdha, M., Miazlovic, H., Foster, T. J., O'Brien, F. J., et al. (2011). Staphylococcus aureus protein A binds to osteoblasts and triggers signals that weaken bone in osteomyelitis. *PLoS One* 6:e18748. doi: 10.1371/journal.pone.0018748

Corrado, A., Donato, P., Maccari, S., Cecchi, R., Spadafina, T., Arcidiacono, L., et al. (2016). Staphylococcus aureus-dependent septic arthritis in murine knee joints: local immune response and beneficial effects of vaccination. *Sci. Rep.* 6:38043. doi: 10.1038/srep38043

Crémieux, A.-C., and Carbon, C. (1997). Experimental models of bone and prosthetic joint infections. *Clin. Infect. Dis.* 25, 1295–1302. doi: 10.1086/516135

Dallas, S. L., and Bonewald, L. F. (2010). Dynamics of the transition from osteoblast to osteocyte. *Ann. N. Y. Acad. Sci.* 1192, 437–443. doi: 10.1111/j.1749-6632.2009.05246.x

Dallas, S. L., Prideaux, M., and Bonewald, L. F. (2013). The osteocyte: an endocrine cell... and more. *Endocr. Rev.* 34, 658–690. doi: 10.1210/er.2012-1026

Deshmane, S. L., Kremlev, S., Amini, S., and Sawaya, B. E. (2009). Monocyte chemoattractant protein-1 (MCP-1): an overview. *J. Interferon Cytokine Res.* 29, 313–326. doi: 10.1089/jir.2008.0027

Eriksen, E. F. (2010). Cellular mechanisms of bone remodeling. *Rev. Endocr. Metab. Disord.* 11, 219–227. doi: 10.1007/s11154-010-9153-1

Florencio-Silva, R., Sasso, G. R., da, S., Sasso-Cerri, E., Simões, M. J., and Cerri, P. S. (2015). Biology of bone tissue: structure, function, and factors that influence bone cells. *Biomed Res. Int.* 2015:421746. doi: 10.1155/2015/421746

Folcik, V. A., An, G. C., and Orosz, C. G. (2007). The Basic Immune Simulator: an agent-based model to study the interactions between innate and adaptive immunity. *Theor. Biol. Med. Model.* 4:39.

Franz-Odendaal, T. A., Hall, B. K., and Witten, P. E. (2006). Buried alive: How osteoblasts become osteocytes. *Dev. Dyn.* 235, 176–190. doi: 10.1002/dvdy.20603

Fullilove, S., Jellis, J., Hughes, S. P. F., Remick, D. G., and Friedland, J. S. (2000). Local and systemic concentrations of tumour necrosis factor-α, interleukin-6 and interleukin-8 in bacterial osteomyelitis. *Trans. R. Soc. Trop. Med. Hyg.* 94, 221–224. doi: 10.1016/S0035-9203(00)90284-0

Fux, C. A., Costerton, J. W., Stewart, P. S., and Stoodley, P. (2005). Survival strategies of infectious biofilms. *Trends Microbiol.* 13, 34–40. doi: 10.1016/j.tim.2004.11.010

Gelfand, A. (2013). The biology of interacting things: the intuitive power of agent-based models. *Biomed. Comput. Rev.* 1, 20–27.

Gentry, L. O. (1997). Management of osteomyelitis. *Int. J. Antimicrob. Agents* 9, 37–42.

Ghosh, S., Matsuoka, Y., Asai, Y., Hsin, K.-Y., and Kitano, H. (2011). Software for systems biology: from tools to integrated platforms. *Nat. Rev. Genet.* 12, 821–832. doi: 10.1038/nrg3096

Ginhoux, F., and Jung, S. (2014). Monocytes and macrophages: developmental pathways and tissue homeostasis. *Nat. Rev. Immunol.* 14, 392–404. doi: 10.1038/nri3671

Goggin, P., Zygalakis, K., Oreffo, R., and Schneider, P. (2016). High-resolution 3D imaging of osteocytes and computational modelling in mechanobiology: insights on bone development, ageing, health and disease. *Eur. Cells Mater.* 31, 264–295. doi: 10.22203/ecm.v031a18

Gong, D., Shi, W., Yi, S., Chen, H., Groffen, J., and Heisterkamp, N. (2012). TGFβ signaling plays a critical role in promoting alternative macrophage activation. *BMC Immunol.* 13:31. doi: 10.1186/1471-2172-13-31

Gopalakrishnan, V., Kim, M., and An, G. (2013). Using an agent-based model to examine the role of dynamic bacterial virulence potential in the pathogenesis of surgical site infection. *Adv. Wound Care* 2, 510–526. doi: 10.1089/wound.2012.0400

Grammatico-Guillon, L., Baron, S., Gettner, S., Lecuyer, A.-I., Gaborit, C., Rosset, P., et al. (2012). Bone and joint infections in hospitalized patients in France, 2008: clinical and economic outcomes. *J. Hosp. Infect.* 82, 40–48. doi: 10.1016/j.jhin.2012.04.025

Hammond, R. A. (2015). *Considerations and Best Practices in Agent-Based Modeling to Inform Policy*. Washington, DC: National Academies Press.

Hernán, M. A. (2015). Invited commentary: agent-based models for causal inference—reweighting data and theory in epidemiology. *Am. J. Epidemiol.* 181, 103–105. doi: 10.1093/aje/kwu272

Holcombe, M., Adra, S., Bicak, M., Chin, S., Coakley, S. I, Graham, A., et al. (2012). Modelling complex biological systems using an agent-based approach. *Integr. Biol.* 4, 53–64. doi: 10.1039/C11B00042J

Horn, F., Heinekamp, T., Kniemeyer, O., Pollmächer, J., Valiante, V., and Brakhage, A. A. (2012). Systems biology of fungal infection. *Front. Microbiol.* 3:108. doi: 10.3389/fmicb.2012.00108

Italiani, P., and Boraschi, D. (2014). From monocytes to M1/M2 macrophages: phenotypical vs. functional differentiation. *Front. Immunol.* 5:514. doi: 10.3389/fimmu.2014.00514

Jilka, R. L., Weinstein, R. S., Bellido, T., Parfitt, A. M., and Manolagas, S. C. (1998). Osteoblast programmed cell death (Apoptosis): modulation by growth factors and cytokines. *J. Bone Miner Res.* 13, 793–802. doi: 10.1359/jbmr.1998.13.5.793

Johansen, L. K., Svalastoga, E. L., Frees, D., Aalbæk, B., Koch, J., Iburg, T. M., et al. (2013). A new technique for modeling of hematogenous osteomyelitis in pigs: inoculation into femoral artery. *J. Invest. Surg.* 26, 149–153. doi: 10.3109/08941939.2012.718043

Jorge, L. S., Chueire, A. G., and Rossit, A. R. B. (2010). Osteomyelitis: a current challenge. *Braz. J. Infect. Dis.* 14, 310–315. doi: 10.1590/S1413-86702010000300020
Josse, J., Velard, F., and Gangloff, S. C. (2015). Staphylococcus aureus vs. Osteoblast: relationship and consequences in osteomyelitis. *Front. Cell. Infect. Microbiol.* 5:85. doi: 10.3389/fcimb.2015.00085

Kaplanski, G., Marin, V., Montero-Julian, F., Mantovani, A., and Farnarier, C. (2003). IL-6: a regulator of the transition from neutrophil to monocyte recruitment during inflammation. *Trends Immunol.* 24, 25–29. doi: 10.1016/S1471-4906(02)00013-3

Kaufmann, S. H. E., and Dorhoi, A. (2016). Molecular determinants in phagocyte-bacteria interactions. *Immunity* 44, 476–491. doi: 10.1016/j.immuni.2016.02.014

Kim, N. D., and Luster, A. D. (2015). The role of tissue resident cells in neutrophil recruitment. *Trends Immunol.* 36, 547–555. doi: 10.1016/j.it.2015.07.007

Kitano, H. (2002a). Computational systems biology. *Nature* 420, 206–210. doi: 10.1038/nature01254

Kitano, H. (2002b). Systems biology: a brief overview. *Science* 295, 1662–1664. doi: 10.1126/science.1069492

Knapp, S., Thalhammer, F., Locker, G. J., Laczika, K., Hollenstein, U., Frass, M., et al. (1998). Prognostic value of MIP-1α, TGF-β2, sELAM-1, and sVCAM-1 in patients with gram-positive sepsis. *Clin. Immunol. Immunopathol.* 87, 139–144. doi: 10.1006/clin.1998.4523

Komarova, S. V., Smith, R. J., Dixon, S. J., Sims, S. M., and Wahl, L. M. (2003). Mathematical model predicts a critical role for osteoclast autocrine regulation in the control of bone remodeling. *Bone* 33, 206–215. doi: 10.1016/S8756-3282(03)00157-1

Kremers, H. M., Nwojo, M. E., Ransom, J. E., Wood-Wentz, C. M., Melton, L. J., and Huddleston, P. M. (2015). Trends in the epidemiology of osteomyelitis: a population-based study, 1969 to 2009. *J. Bone Joint Surg. Am.* 97, 837–845. doi: 10.2106/JBJS.N.01350

Kurtz, S. M., Ong, K. L., Lau, E., Widmer, M., Maravic, M., Gómez-Barrena, E., et al. (2011). International survey of primary and revision total knee replacement. *Int. Orthop.* 35, 1783–1789. doi: 10.1007/s00264-011-1235-5

Laurent, E., Gras, G., Druon, J., Rosset, P., Bernard, L., and Grammatico-Guillon, L. (2016). Évolution des infections ostéo-articulaires (IOA) en France après mise en place des centres de références des IOA complexes (CRIOAC) : PMSI 2008 versus 2013. *Rev. d'Épidémiol. Santé Publ.* 64, S23–S24. doi: 10.1016/j.respe.2016.01.075

Lebeaux, D., Chauhan, A., Rendueles, O., and Beloin, C. (2013). From *in vitro* to *in vivo* models of bacterial biofilm-related infections. *Pathogens* 2, 288–356. doi: 10.3390/pathogens2020288

Lee, Y.-W., Chung, Y., Juhn, S. K., Kim, Y., and Lin, J. (2011). Activation of the transforming growth factor-beta pathway in bacterial otitis media. *Ann. Otol. Rhinol. Laryngol.* 120, 204–213. doi: 10.1177/00034894111200310

Lio, P., Paoletti, N., Moni, M. A., Atwell, K., Merelli, E., and Viceconti, M. (2012). Modelling osteomyelitis. *BMC Bioinformatics* 13(Suppl. 14):S12. doi: 10.1186/1471-2105-13-S14-S12

Lipsky, B., Weigelt, J., Gupta, V., Killian, A., and Peng, M. (2007). Skin, soft tissue, bone, and joint infections in hospitalized patients: epidemiology and microbiological, clinical, and economic outcomes. *Infect. Control Hosp. Epidemiol.* 28, 1290–1298. doi: 10.1086/520743

Mader, J. T. (1985). Animal models of osteomyelitis. *Am. J. Med.* 78, 213–217. doi: 10.1016/0002-9343(85)90387-0

Manolagas, S. C. (2000). Birth and death of bone cells: basic regulatory mechanisms and implications for the pathogenesis and treatment of osteoporosis. *Endocr. Rev.* 21, 115–137. doi: 10.1210/edrv.21.2.0395

Marriott, I., Gray, D. L., Rati, D. M., Fowler, V. G., Stryjewski, M. E., Levin, L. S., et al. (2005). Osteoblasts produce monocyte chemoattractant protein-1 in a murine model of Staphylococcus aureus osteomyelitis and infected human bone tissue. *Bone* 37, 504–512. doi: 10.1016/j.bone.2005.05.011

Olson, M. E., and Horswill, A. R. (2013). *Staphylococcus aureus* Osteomyelitis: bad to the bone. *Cell Host Microbe* 13, 629–631. doi: 10.1016/j.chom.2013.05.015

Paoletti, N., Lio, P., Merelli, E., and Viceconti, M. (2012). Multilevel computational modeling and quantitative analysis of bone remodeling.

IEEE/ACM Trans. Comput. Biol. Bioinform. 9, 1366–1378. doi: 10.1109/TCBB.2012.51

Parfitt, A. M. (1994). Osteonal and hemi-osteonal remodeling: the spatial and temporal framework for signal traffic in adult human bone. *J. Cell. Biochem.* 55, 273–286. doi: 10.1002/jcb.240550303

Parwaresch, M. R., and Wacker, H. H. (1984). Origin and kinetics of resident tissue macrophages. Parabiosis studies with radiolabelled leucocytes. *Cell Tissue Kinet.* 17, 25–39. doi: 10.1111/j.1365-2184.1984.tb00565.x

Patel, A. A., Zhang, Y., Fullerton, J. N., Boelen, L., Rongvaux, A., Maini, A. A., et al. (2017). The fate and lifespan of human monocyte subsets in steady state and systemic inflammation. *J. Exp. Med.* 214, 1913–1923. doi: 10.1084/jem.20170355

Politopoulos, I. (2007). Review and analysis of agent-based models in biology. *Proc. Natl. Acad. Sci. U.S.A.* 104, 12017–12022.

Prame Kumar, K., Nicholls, A. J., and Wong, C. H. Y. (2018). Partners in crime: neutrophils and monocytes/macrophages in inflammation and disease. *Cell Tissue Res.* 371, 551–565. doi: 10.1007/s00441-017-2753-2

Rankin, S. M. (2010). The bone marrow: a site of neutrophil clearance. *J. Leukoc. Biol.* 88, 241–251. doi: 10.1189/jlb.0210112

Repp, F., Kollmannsberger, P., Roschger, A., Kerschnitzki, M., Berzlanovich, A., Gruber, G. M., et al. (2017). Spatial heterogeneity in the canalicular density of the osteocyte network in human osteons. *Bone Rep.* 6, 101–108. doi: 10.1016/j.bonr.2017.03.001

Roodman, G. D. (1999). Cell biology of the osteoclast. *Exp. Hematol.* 27, 1229–1241. doi: 10.1016/S0301-472X(99)00061-2

Ryser, M. D., Nigam, N., and Komarova, S. V. (2009). Mathematical modeling of spatio-temporal dynamics of a single bone multicellular unit. *J. Bone Miner. Res.* 24, 860–870. doi: 10.1359/jbmr.081229

Schneider, P., Meier, M., Wepf, R., and Müller, R. (2010). Towards quantitative 3D imaging of the osteocyte lacuno-canalicular network. *Bone* 47, 848–858. doi: 10.1016/j.bone.2010.07.026

Shapiro, F. (2008). Bone development and its relation to fracture repair. The role of mesenchymal osteoblasts and surface osteoblasts. *Eur. Cell Mater.* 15, 53–76. doi: 10.22203/ecm.v015a05

Shi, Z. Z., Wu, C.-H., and Ben-Arieh, D. (2014). Agent-based model: a surging tool to simulate infectious diseases in the immune system. *Open J. Model. Simul.* 02, 12–22. doi: 10.4236/ocmsi.2014.21004

Smith, A. M., McCullers, J. A., and Adler, F. R. (2011). Mathematical model of a three-stage innate immune response to a pneumococcal lung infection. *J. Theor. Biol.* 276, 106–116. doi: 10.1016/j.jtbi.2011.01.052

Smith, I. M., Austin, O. M. B., and Batchelor, A. G. (2006). The treatment of chronic osteomyelitis: a 10 year audit. *J. Plast. Reconstr. Aesthet. Surg.* 59, 11–15. doi: 10.1016/j.bjps.2005.07.002

Stengel, D., Bauwens, K., Sehouli, J., Ekkernkamp, A., and Porzolt, F. (2001). Systematic review and meta-analysis of antibiotic therapy for bone and joint infections. *Lancet Infect. Dis.* 1, 175–188. doi: 10.1016/S1473-3099(01)00094-9

Strydom, N., and Rankin, S. M. (2013). Regulation of circulating neutrophil numbers under homeostasis and in disease. *J. Innate Immun.* 5, 304–314. doi: 10.1159/000350282

Tecchio, C., Micheletti, A., and Cassatella, M. A. (2014). Neutrophil-derived cytokines: facts beyond expression. *Front. Immunol.* 5:08. doi: 10.3389/fimmu.2014.00508

Vashishth, D., Verborgt, O., Divine, G., Schaffler, M. B., and Fyhrie, D. P. (2000). Decline in osteocyte lacunar density in human cortical bone is associated with accumulation of microcracks with age. *Bone* 26, 375–380. doi: 10.1016/S8756-3282(00)00236-2

Wagner, C., Kondella, K., Bernschneider, T., Heppert, V., Wentzensen, A., and Hänsch, G. M. (2003). Post-traumatic osteomyelitis: analysis of inflammatory cells recruited into the site of infection. *Shock* 20, 503–510. doi: 10.1097/01.shk.0000093542.78705.e3

Wallace, W. A., Gillooly, M., and Lamb, D. (1993). Age related increase in the intra-alveolar macrophage population of non-smokers. *Thorax* 48, 668–669. doi: 10.1136/thx.48.6.668
Walter, G., Kemmerer, M., Kappler, C., and Hoffmann, R. (2012). Treatment algorithms for chronic osteomyelitis. *Dtsch. Arztebl. Int.* 109, 257–264. doi: 10.3238/arztebl.2012.0257

Whitelaw, D. M., and Bell, M. (1966). The intravascular lifespan of monocytes. *Blood* 28, 455–464. doi: 10.1182/blood.v28.3.455.455

Wilensky, U. (1999). *NetLogo. Center for Connected Learning and Computer-Based Modeling*. Available online at: https://ccl.northwestern.edu/netlogo/index.shtml
(accessed October 13, 2017).

Xing, J., Hou, T., Luobu, B., Luo, F., Chen, Q., Li, Z., et al. (2012). Anti-infection tissue engineering construct treating osteomyelitis in rabbit tibia. *Tissue Eng. Part A* 19, 255–263. doi: 10.1089/ten.tea.2012.0262

**Conflict of Interest:** SAI and VL were employed by company Tekliko SARL.

The remaining 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.

Copyright © 2020 Alsassa, Lefèvre, Laugier, Stindel and Ansart. This is an open-
access article distributed under the terms of the Creative Commons Attribution
License (CC BY). The use, distribution or reproduction in other forums is permitted,
provided the original author(s) and the copyright owner(s) are credited and that the
original publication in this journal is cited, in accordance with accepted academic
practice. No use, distribution or reproduction is permitted which does not comply
with these terms.