(C) PLOS One
This story was originally published by PLOS One and is unaltered.
. . . . . . . . . .



Integrated computational and in vivo models reveal Key Insights into macrophage behavior during bone healing [1]

['Etienne Baratchart', 'Integrated Mathematical Oncology Department', 'Moffitt Cancer Center', 'Research Institute', 'Tampa', 'Florida', 'United States Of America', 'Chen Hao Lo', 'Cancer Biology Ph.D. Program', 'Department Of Cell Biology Microbiology']

Date: 2022-07

Data from in vivo bone injury experiment were derived from our previous work [ 40 ]. In this study, all animal studies were performed in accordance with Guidelines for the Care and Use of Laboratory Animals published by the National Institutes of Health, and approved by the Animal Care and Use Committee at the University of South Florida, under IACUC Protocol R5857 (CCL). Male C57BL/6 mice (5–6 weeks old) were purchased from Jackson Laboratory. Mice (n = 30) were subject to tibial bone injury by penetration of a 28-gauge (0.3062mm diameter) syringe through the knee epiphysis to mid-shaft.

Data from in vivo bone injury experiment were derived from our previous work [ 40 ]. In this study, all animal studies were performed in accordance with Guidelines for the Care and Use of Laboratory Animals published by the National Institutes of Health, and approved by the Animal Care and Use Committee at the University of South Florida, under IACUC Protocol R5857 (CCL). Male C57BL/6 mice (5–6 weeks old) were purchased from Jackson Laboratory. Mice (n = 30) were subject to tibial bone injury by penetration of a 28-gauge (0.3062mm diameter) syringe through the knee epiphysis to mid-shaft. Mice tibias at baseline and at days 1, 2, 3, 7 and 14 (n = 5/time point) were collected for analysis. Temporal population data was used to parameterize subsequent mathematical models.

Bone volume data was derived from formalin-fixed tibias by micro-computed topography (μCT) scanning using Scanco μ35 scanner. Endosteal trabecular bone volume was analyzed 100μm away from the tip of growth plate to clear the dense bone nature of the growth plate. 1000μm along the midshaft of each bone was then scanned and analyzed using built-in functions (n = 30 bones; 5/time point).

After μCT analysis, tibia bones were decalcified using 14% EDTA for 3 weeks for further staining quantitation and analyses. Decalcified bones were sectioned at 4μm thickness. Sections were enzymatically stained for tartrate-resistant acid phosphatase (TRAcP) for osteoclast numbers based on manufacturer’s protocol [ 120 ]. Stained slides were imaged using the Evos Auto microscope to capture 20X photos which included injury site and its immediate periphery. All TRAcP positive (red) multinucleated osteoclasts within 5μm radius from injury were counted, and mathematically converted to osteoclasts / bone marrow volume (#OCL/μm 3 ) for each bone at each time point.

FFPE tibia bones were further sequentially sectioned and baked at 56°C in preparation for immunofluorescence staining of osteoblast (RUNX2 at 1:500; Abcam Cat. No. ab81357) and nuclear staining (DAPI). Deparaffined and rehydrated slides were subject to heat-induced antigen retrieval method. Sections were then blocked and incubated in primary antibodies diluted in 10% normal goat serum in TBS overnight at 4°C. Subsequently, slides were stained with secondary Alexa Fluor 568-conjugated antibody at 1:1000 at room temperature for 1 hour under light-proof conditions. Stained slides were stained with DAPI for nuclear contrast and mounted for imaging at 20X using Zeiss upright fluorescent microscope to include the injury site as well as the immediate peripheral tissue. All runx2 positive cells (red staining colocalizing with DAPI) within 5μm radius from injury were counted and mathematically converted to osteoblasts / bone marrow volume (#OBL/μm 3 ) for each bone at each time point.

Harvested contralateral injured tibias (n = 30; 5/time point) had tips removed and were subjected to centrifugation at 16,000g for 5 seconds for isolation of whole bone marrow for flow cytometry staining and analysis. Red blood cells were lysed using RBC Lysis Buffer from Sigma Aldrich (Cat. No. R7757-100ML) as per manufacturer’s guidelines. Live bone marrow cells were subject to FcR-receptor blocking (1:3; BioLegend; Cat. No. 101319) and viability staining (1:500; BioLegend; Cat. No. 423105). Samples were then stained by cell-surface conjugated antibodies from BioLegend diluted in autoMACS buffer (Miltenyi; Cat. No. 130-091-221) for phenotyping myeloid cells: CD11b-BV786 (1:200; Cat. No. 101243), LY-6C-Alexa Fluor 488 (1:500; Cat. No.128021) and LY-6G-Alexa Fluor 700 (1:200; Cat. No. 561236). Cells were then fixed with 2% paraformaldehyde in dark prior to intracellular staining. Fixed cells were permeabilized using intracellular conjugated antibodies to assess polarization status: NOS2-APC (1:100; eBioscience; Cat. No. 17-5920-80) and ARG1-PE (1:100; R&D; Cat. No. IC5868P). Appropriate compensation and fluorescence-minus-one (FMO) controls were generated in parallel either with aliquots of bone marrow cells or Rainbow Fluorescent Particle beads (BD Biosciences; Cat. No. 556291). All antibody concentrations were titrated prior to injury study using primary bone marrow cells to ensure optimal separation and detection of true negative and positive populations. Stained controls and samples were analyzed using BD Biosciences LSR flow cytometer ( S2 Fig ).

Mathematical and computational methods

Comprehensive model structure. In order to efficiently describe the process of model building regarding the different hypothesis combinations (Fig 2), we begin with a generic set of coupled equations, using a formulation that is valid for all hypotheses combinations (S3 Fig). For all populations, homeostasis was described by a replenishment term and a clearance term. For polarized monocytes/ macrophages, no replenishment was considered at homeostasis, as the baseline measured by flow cytometry was close to zero. Here is the detailed description, equation by equation:

Equation 1: Naive monocytes. Naive monocytes are assumed to be replenished at a constant rate H Mo , and to die at a rate δ Mo Mo. δ Mo , the lifespan parameter, was retrieved from literature (Table 4), and H Mo was estimated so that monocyte level at homeostasis would match experimentally measured monocyte baseline. The term I(M 1 , Mo 1 , D) corresponds to the number of monocytes infiltrating the bone marrow per unit of time due to injury factors and pro-inflammatory cells. It is equal to I 1 (M 1 + Mo 1 ) + I 2 D. As indicated by the mathematical formulation, inflammation-associated monocyte recruitment is driven by injury signals on one hand, and pro-inflammatory monocytes/macrophages on the other hand. This reflects the fact that cellular debris and pro-inflammatory cells produce factors that help recruiting additional monocytes [1]. The term p 3 (D, M 1 , Mo 1 )Mo represents the number of naive monocytes polarizing into pro-inflammatory monocytes per unit of time, as a function of cellular debris and pro-inflammatory cells. The term depol 3 (M 2 )Mo 1 represents the number pro-inflammatory monocytes reverting back to a naive state per unit of time, as function of anti-inflammatory macrophages.

Equation 2: Naive macrophages. Naive macrophage are assumed to be replenished at a constant rate H M , and to die at a rate δ M M. δ M , the lifespan parameter, was retrieved from literature (Table 4), and H M was estimated so that macrophage level at homeostasis would match experimentally measured macrophage baseline. The term d OC (OB, M 1 , Mo 1 , M 2 )M corresponds to the number of macrophages differentiating into osteoclasts per unit of time, as a function of osteoblasts (RANKL, OPG), pro-inflammatory macrophages (IL-1, TNF), pro-inflammatory monocytes (IL-1, TNF) and anti-inflammatory macrophages (IL-10, TGF). The term p 1 (D, M 1 , Mo 1 )M represents the number of naive macrophages polarizing into a pro-inflammatory state per unit of time, as a function of injury factors and pro-inflammatory cells. The term p 21 (D, M 1 , Mo 1 )M represents the number of naive macrophages polarizing into an anti-inflammatory state per unit of time, as a function of injury factors and pro-inflammatory cells. The term depol 1 (M 2 )M 1 represents the number of pro-inflammatory macrophages that revert back to a naive state per unit of time, as a function of anti-inflammatory macrophages. The term depol 21 M 2 represents the number of anti-inflammatory macrophages that revert back to a naive state per unit of time. Of note, no influx or differentiation from monocytes during inflammation were assumed, as macrophage population does not show evidence of expansion in our obtained biological data.

Equation 3: Pro-inflammatory macrophages. Pro-inflammatory macrophages are assumed to be absent of the bone marrow under homeostatic conditions, as experimental al baseline level did not 0.3% across the replicates. The term p 1 (D, M 1 , Mo 1 )M represents the number of pro-inflammatory macrophages generated (from naive pool) per unit of time, as a function of injury factors and pro-inflammatory cells. The term depol 1 (M 2 )M 1 represents the number of pro-inflammatory macrophages that revert back to a naive state per unit of time, as a function of anti-inflammatory macrophages. The term depol 22 M 2 represents the number of anti-inflammatory macrophages that reprogram into a pro-inflammatory phenotype. The term p 22 (D)M 1 represents the number of pro-inflammatory macrophages that reprogram into an anti-inflammatory phenotype, as a function of injury factors.

Equation 4: Anti-inflammatory macrophages. Anti-inflammatory macrophages are assumed to be absent of the bone marrow under homeostatic conditions, as experimental al baseline level did not exceed 0.07% across the replicates. The term p 21 (D, M 1 , Mo 1 )M represents the number of anti-inflammatory macrophages generated (from naive pool) per unit of time, as a function of injury factors and pro-inflammatory cells. The term p 22 (D)M 1 represents the number of pro-inflammatory macrophages that reprogram into an anti-inflammatory phenotype, as a function of injury factors. The term depol 21 M 2 represents the number of pro-inflammatory macrophages that revert back to a naive state per unit of time. The term depol 22 M 2 represents the number of anti-inflammatory macrophages that reprogram into a pro-inflammatory phenotype.

Equation 5: Pro-inflammatory monocytes. Pro-inflammatory monocytes are assumed to be absent of the bone marrow under homeostatic conditions, as experimental al baseline level did not exceed 0.16% across the replicates. The term p 3 (D, M 1 , Mo 1 )Mo represents the number of pro-inflammatory monocytes generated (from naive pool) per unit of time, as a function of injury factors and pro-inflammatory cells. The term depol 3 (M 2 )Mo 1 represents the number of pro-inflammatory monocytes that revert back to a naive state per unit of time, as a function of anti-inflammatory macrophages. The term depol 22 M 2 represents the number of anti-inflammatory macrophages that reprogram into a pro-inflammatory phenotype. The term p 22 (D)M 1 represents the number of pro-inflammatory macrophages that reprogram into an anti-inflammatory phenotype, as a function of injury factors.

Pro-inflammatory macrophages/monocytes and anti-inflammatory macrophages: Hypothesis c1. Both CD11b+Ly6C+ monocytes and CD11b+Ly6C- macrophages can polarize into a pro-inflammatory phenotype [1], in response to injury signals or to factors produced by already present pro-inflammatory cells [1]. In assumptions c1 and c2, monocytes and macrophages polarize into pro-inflammatory monocytes and macrophages respectively through two terms. The first is proportional to the amount of injury signals present and the second is proportional to the amount of pro-inflammatory cells present. In the assumption c3, injury signals polarize local resident macrophages, and those in turn promote polarization of pro-inflammatory monocytes. In c3, pro-inflammatory macrophages repolarize into anti-inflammatory macrophages by the uptake of cellular debris/injury signals [121]. In this scenario, by plasticity, anti-inflammatory macrophages naturally return to a pro-inflammatory phenotype in absence of signals [122]. Polarization into anti-inflammatory macrophages was assumed to be proportional to the amount of injury signals for assumption c1, proportional to the amount of pro-inflammatory cells for assumption c2, and a transition term from pro-inflammatory macrophages for c3.

Pro-inflammatory macrophages/monocytes and anti-inflammatory macrophages: Hypothesis c2. Both CD11b+Ly6C+ monocytes and CD11b+Ly6C- macrophages can polarize into a pro-inflammatory phenotype [1], in response to injury signals or to factors produced by already present pro-inflammatory cells [1]. In assumptions c1 and c2, monocytes and macrophages polarize into pro-inflammatory monocytes and macrophages respectively through two terms. The first is proportional to the amount of injury signals present and the second is proportional to the amount of pro-inflammatory cells present. In the assumption c3, injury signals polarize local resident macrophages, and those in turn promote polarization of pro-inflammatory monocytes. In c3, pro-inflammatory macrophages repolarize into anti-inflammatory macrophages by the uptake of cellular debris/injury signals [121]. In this scenario, by plasticity, anti-inflammatory macrophages naturally return to a pro-inflammatory phenotype in absence of signals [122]. Polarization into anti-inflammatory macrophages was assumed to be proportional to the amount of injury signals for assumption c1, proportional to the amount of pro-inflammatory cells for assumption c2, and a transition term from pro-inflammatory macrophages for c3.

Pro-inflammatory Macrophages/Monocytes and Anti-inflammatory Macrophages: Hypothesis c3. Both CD11b+Ly6C+ monocytes and CD11b+Ly6C- macrophages can polarize into a pro-inflammatory phenotype [1], in response to injury signals or to factors produced by already present pro-inflammatory cells [1]. In assumptions c1 and c2, monocytes and macrophages polarize into pro-inflammatory monocytes and macrophages respectively through two terms. The first is proportional to the amount of injury signals present and the second is proportional to the amount of pro-inflammatory cells present. In the assumption c3, injury signals polarize local resident macrophages, and those in turn promote polarization of pro-inflammatory monocytes. In c3, pro-inflammatory macrophages repolarize into anti-inflammatory macrophages by the uptake of cellular debris/injury signals [121]. In this scenario, by plasticity, anti-inflammatory macrophages naturally return to a pro-inflammatory phenotype in absence of signals [122]. Polarization into anti-inflammatory macrophages was assumed to be proportional to the amount of injury signals for assumption c1, proportional to the amount of pro-inflammatory cells for assumption c2, and a transition term from pro-inflammatory macrophages for c3.

Equation 6: Osteoblasts. Osteoblast are assumed to be replenished at a rate H OB OC, proportional to osteoclasts. This reflects the ability of osteoclasts to produce osteogenic signals like transforming growth factor β (TGFβ) and bone morphogenetic proteins (BMPs) [123]. Similar assumptions are considered in published works of homeostatic bone remodeling [32–34]. Osteoblasts are assumed to die at a rate δ OB OBB. δ OB , the lifespan parameter, was retrieved from literature (Table 4), and H OB was estimated so that osteoblast level at homeostasis would match experimentally measured osteoblast baseline. The term γ OB (D, M 2 ) represents the number of osteoblasts generated per unit of time, as a function of injury factors and anti-inflammatory factors.

Osteoblast dynamics: Hypothesis b1. The osteoblast clearance term was assumed to be proportional to the bone volume, in order to account for osteoblast differentiation into osteocytes, when resorbing bone matrix [123]. A similar assumption is made in the model developed by Ryser et al. that describes bone remodeling as a spatial evolutionary game [124]. During injury, an extra term for osteoblast expansion is present, driven by anti-inflammatory macrophages (hypothesis b1) or injury factors (hypothesis b2), both supported by literature [1,3,6,7,125].

Osteoblast dynamics: Hypothesis b2. The osteoblast clearance term was assumed to be proportional to the bone volume, in order to account for osteoblast differentiation into osteocytes, when resorbing bone matrix [123]. A similar assumption is made in the model developed by Ryser et al. that describes bone remodeling as a spatial evolutionary game [124]. During injury, an extra term for osteoblast expansion is present, driven by anti-inflammatory macrophages (hypothesis b1) or injury factors (hypothesis b2), both supported by literature [1,3,6,7,125].

Equation 7: Osteoclasts. Osteoclasts are assumed to be replenished at a rate d OC (OB, M 1 , Mo 1 , M 2 )M, which reflects differentiation of macrophages into osteoclasts, as a function of osteoblasts, pro-inflammatory macrophages, pro-inflammatory monocytes, anti-inflammatory macrophages [126–128]. This reflects osteoclastic factors produced by osteoblasts (RANKL) and pro-inflammatory monocytes/macrophages (IL-1, TNF), as well as anti-osteoclastic factors produced by osteoblasts (OPG) and anti-inflammatory macrophages (transforming growth factor β (TGFβ), IL-10). The term δ OC (M 2 , OB)OC represents the number of osteoclasts dying per unit of time, as a function of anti-inflammatory macrophages and osteoblasts. This reflects factors produced by anti-inflammatory macrophages (IL-10, TGFβ) and osteoblasts (OPG) that reduce osteoclast lifespan.

Osteoclast dynamics: Hypothesis a1. This osteoclast formation term was assumed to be proportional to osteoblasts for assumption a3, reflecting the ability of osteoblastic cells to produce RANKL, which is an essential mediator of osteoclast formation [123]. For the other assumptions, homeostatic osteoclast replenishment was assumed to be constant. This term had an additional contribution from pro-inflammatory monocytes/macrophages for assumptions a1 and a2, representing the ability of pro-inflammatory monocytes/macrophages to produce factors like IL-1 and TNF that favor osteoclast formation [123,129]. Osteoclast formation was divided by an inhibitory term, a linear function of anti-inflammatory macrophages for assumptions a1 and a3, and a linear function of osteoblasts for assumption a2. The first assumption reflects factors produced by anti-inflammatory macrophages, like IL-10, that disrupt osteoclast formation [86], whereas the second reflects the ability of osteoblasts to produce osteoprotegerin (OPG), a RANKL decoy receptor [123]. Moreover, this inhibition affects not only the ability of monocytes-macrophages to fuse and form osteoclasts, but also their life span. Indeed RANKL is necessary for osteoclast survival since OPG produced by osteoblasts reduces their life span [130]. Similarly, anti-inflammatory macrophages produce TGFβ, which is known to drive osteoclast apoptosis [131].

Osteoclast dynamics: Hypothesis a2. This osteoclast formation term was assumed to be proportional to osteoblasts for assumption a3, reflecting the ability of osteoblastic cells to produce RANKL, which is an essential mediator of osteoclast formation [123]. For the other assumptions, homeostatic osteoclast replenishment was assumed to be constant. This term had an additional contribution from pro-inflammatory monocytes/macrophages for assumptions a1 and a2, representing the ability of pro-inflammatory monocytes/macrophages to produce factors like IL-1 and TNF that favor osteoclast formation [123,129]. Osteoclast formation was divided by an inhibitory term, a linear function of anti-inflammatory macrophages for assumptions a1 and a3, and a linear function of osteoblasts for assumption a2. The first assumption reflects factors produced by anti-inflammatory macrophages, like IL-10, that disrupt osteoclast formation [86], whereas the second reflects the ability of osteoblasts to produce osteoprotegerin (OPG), a RANKL decoy receptor [123]. Moreover, this inhibition affects not only the ability of monocytes-macrophages to fuse and form osteoclasts, but also their life span. Indeed RANKL is necessary for osteoclast survival since OPG produced by osteoblasts reduces their life span [130]. Similarly, anti-inflammatory macrophages produce TGFβ, which is known to drive osteoclast apoptosis [131].

Osteoclast dynamics: Hypothesis a3. This osteoclast formation term was assumed to be proportional to osteoblasts for assumption a3, reflecting the ability of osteoblastic cells to produce RANKL, which is an essential mediator of osteoclast formation [123]. For the other assumptions, homeostatic osteoclast replenishment was assumed to be constant. This term had an additional contribution from pro-inflammatory monocytes/macrophages for assumptions a1 and a2, representing the ability of pro-inflammatory monocytes/macrophages to produce factors like IL-1 and TNF that favor osteoclast formation [123,129]. Osteoclast formation was divided by an inhibitory term, a linear function of anti-inflammatory macrophages for assumptions a1 and a3, and a linear function of osteoblasts for assumption a2. The first assumption reflects factors produced by anti-inflammatory macrophages, like IL-10, that disrupt osteoclast formation [86], whereas the second reflects the ability of osteoblasts to produce osteoprotegerin (OPG), a RANKL decoy receptor [123]. Moreover, this inhibition affects not only the ability of monocytes-macrophages to fuse and form osteoclasts, but also their life span. Indeed RANKL is necessary for osteoclast survival since OPG produced by osteoblasts reduces their life span [130]. Similarly, anti-inflammatory macrophages produce TGFβ, which is known to drive osteoclast apoptosis [131].

Equation 8: Bone volume. Bone dynamics is described by two terms: a bone resorption term, which is the volume of bone resorbed per unit of time and is assumed to be proportional to the number of osteoclasts, and a bone formation term, which is the volume of bone formed per unit of time and is assumed to be proportional to the number of osteoblasts. Such assumptions have been broadly used across a large variety of modeling studies [32–34]. δ B (1 + α)OCB is the bone resorption term and is the sum of the two terms δ B OCB (homeostatic resorption) and α(M 1 + Mo 1 )δ B OCB (pro-inflammatory monocyte/macrophage mediated resorption). As indicated by this mathematical formulation, bone resorption was assumed to also be proportional to the bone mass. This reflects the fact that more bone volume increases the likelihood for bone resorption. Furthermore, this formulation ensures bone mass stays strictly positive in the model. Π B (1 + β)OB is the bone formation term and is the sum of the two terms Π B OB (homeostatic bone formation) and βM 2 Π B OB (anti-inflammatory macrophage mediated bone formation). Under homeostasis, bone is formed at rate Π B OB and resorbed at rate δ B OCB. Resorption rate δ B , as well as as parameters α and β, were calibrated on bone volume dynamics, and bone formation rate Π B was then imposed by the relation Π B OB 0 = δ B OC 0 B 0 , which ensures that bone volume remains at homeostasis when osteoblasts and osteoclasts are at homeostatic levels. As indicated by mathematical formulations, bone formation and resorption terms were assumed to linearly increase with respect to anti-inflammatory and pro-inflammatory monocytes/ macrophages, respectively. This accounts for the fact that anti-inflammatory macrophages typically produce osteogenic factors like TGFβ or OSM, that are known to promote osteoblast expansion and bone mineralization [3], and for the fact that pro-inflammatory monocytes/macrophages typically produce osteolytic factors like TNF and IL-1, that are known to promote osteoclast resorptive activity [132].

Equation 9: Injury factors. Injury factors dynamics consists in an exponential type of decay, with a decay rate δ D (M 1 + Mo 1 ) proportional to pro-inflammatory monocytes/macrophages number, which represents how pro-inflammatory monocytes/macrophages uptake cellular debris, which in return reduces pro-inflammatory signals.

ODE solver. The ODE45 function of Matlab was used to solve the differential equation system. The experimental baseline values (time 0) were used as initial conditions.

Parameter estimation method. To estimate parameters for goodness of fit, we defined the following objective function: Where i represents the time point index and j the variable index, α represents the parameter set used to evaluate the model function f, D ij represents the experimental data of variable j at time point i, σ i represents the experimental standard deviation (over all the animals of a given time point), and N represents the number of time points. The functional J 2 corresponds to the weighted least squares criterion. The functional J ∞ is the Tchebychev approximation, which considers the maximal residual instead of the sum of the residuals. In both cases, the choice of the max over the observed variables of the sum of the squares of the residuals was motivated to make sure that every variable was fitted with equal relative importance. Indeed, in the case of the minimization of the sum of the squares over all the variables, it is sometimes possible to find an optimum in minimizing certain variables at the detriment of others. This way, we ensure that all variables are equally well fitted. The reason for considering the objective function J ∞ in addition to the classical criterion J 2 is to avoid neglecting any time point in the fit. The least squares functional allows sometimes to find an optimum optimizing certain time points at the detriment of others. In this current study, this is a potentially big issue, as the time sampling is not homogenous across the time points, meaning that biological dynamics of importance might be ignored, while still producing a good fit under the least squares metric. In order to minimize this function representing the error estimate between data and model, we used the Matlab function fminsearch with a penalization term to stay in a parameter range set with reasonable boundaries. In order to rank the models in term of goodness of fit, we used AIC, that is defined as follows: where p is the number of parameters, and the functional J is either J 2 or J ∞ .

[END]
---
[1] Url: https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1009839

Published and (C) by PLOS One
Content appears here under this condition or license: Creative Commons - Attribution BY 4.0.

via Magical.Fish Gopher News Feeds:
gopher://magical.fish/1/feeds/news/plosone/

You are viewing proxied material from magical.fish. The copyright of proxied material belongs to its original authors. Any comments or complaints in relation to proxied material should be directed to the original authors of the content concerned. Please see the disclaimer for more details.