Quantitative assessment of extracellular IL-1 regulation

Journal of Theoretical Biology 487 (2019) 110113

Contents lists available at ScienceDirect

Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/jtb

Quantitative assessment of extracellular IL-1 regulation Luis F. Ponce∗, Karina García-Martínez, Kalet León System Biology Department, Center of Molecular Immunology, Habana 11600, Cuba

a r t i c l e

i n f o

Article history: Received 7 May 2019 Revised 10 November 2019 Accepted 8 December 2019 Available online 9 December 2019 Keywords: Interleukin 1 Signaling regulation Mathematical model Soluble interleukin 1 receptors Cell membrane interactions

a b s t r a c t IL-1 system is involved in the induction and maintenance of chronic inflammation associated with several autoimmune diseases and cancer, mainly due to its capacity to promote the secretion of inflammatory mediators. For this reason, several intracellular and extracellular mechanisms for this system have been fixed during the evolution. In spite of the large description of molecular interactions between IL-1 ligands and receptors, little is known about the relevance and limits of the extracellular regulatory mechanims in different scenarios. To tackle this problem, we developed and calibrated a mathematical model including all the known interactions between IL-1 ligands and IL-1Rs and calibrate it with experimental data of IL-1 binding to different cells. The model predicts that, independently on the IL-1Rs expression, IL-1α has more ability than IL-1β to induce IL-1 signaling, which suggests that both ligands can be equally relevant for the IL-1 related inflammation. On the other hand, at the cell level, IL-1 signaling is mainly controlled by IL-1R1 and IL-1R3 and not by IL-1R2. Moreover, the soluble form of IL-1R1 and IL-1RA have the highest capacity to prevent IL-1α while IL-1R2 and IL-1R1 and IL-1RA have a similar capacity to prevent IL-1β signaling. The soluble IL-1R3 has the lowest capacity to prevent IL-1 signaling and preferentially inhibits cells with low number of IL-1R3. In general, model predictions suggest several ways in which IL-1 controlling system may fail, developing IL-1 related inflammation. © 2019 Elsevier Ltd. All rights reserved.

1. Introduction Interleukin 1 (IL-1) is a system of cytokines and receptors involved in the activation of innate and adaptive immune system primarily through the induction of inflammatory mediators (Dinarello, 2009). It plays a major role in the development of several autoimmune diseases such as gout (Pope and Tschopp, 2007), type 2 diabetes (Banerjee and Saxena, 2012), rheumatoid arthritis (Kay and Calaberase, 2004) and Alzheimer (Shaftel et al., 2008). On the other hand, IL-1 system is also involved in the development and invasiveness of several lethal cancers such as lung (Saijo et al., 2002), esophageal (Chen et al., 2012), pancreatic (Zhuang et al., 2016), cervical (Ricote et al., 2004) and colorectal cancer (West et al., 2015). It has been demonstrated that IL-1 promotes tumorigenesis through the induction of processes such as angiogenesis (Voronov et al., 2014), the recruitment of myeloid derived suppressor cells (Tu et al., 2008), and macrophages differentiation to M2 phenotype (Hagemann et al., 2008). For these reasons, blocking IL-1 has been proven to be effective as a therapeutic approach for both cancer and autoimmune diseases (reviewed in (Dinarello, 2010; Dinarello and van der Meer, 2013)).

∗

Corresponding author. E-mail address: [email protected] (L.F. Ponce).

https://doi.org/10.1016/j.jtbi.2019.110113 0022-5193/© 2019 Elsevier Ltd. All rights reserved.

IL-1 system includes three soluble ligands: two agonists known as IL-1α and IL-1β , and a natural antagonist called IL-1 receptor antagonist (IL-1RA). From solution, the agonist ligands bind to two membrane receptors: IL-1 receptor type I (IL-1R1) and IL-1 receptor type II (IL-1R2 also called the decoy receptor) (Garlanda et al., 2013). The resulting dimeric complex from these interactions binds to the IL-1R accessory chain protein (IL-1R3), forming a trimeric complex. Only the complexes formed by IL-1α or IL-1β , IL-1R1 and IL-1R3 are capable to activate the IL-1 singnaling, been responsible for its biological effect (Cullinan et al., 1998). The IL-1RA binds preferentially to IL-1R1, but does not mediate the interaction with IL-1R3. By competing for the IL-1R1, IL-1RA inhibits IL-1 signal and controls the inflammatory processes at the tissue level in several localizations (Arend and Gabay 20 0 0). In addition to the number of IL-1 system components, other aspects can make more complex the regulation of IL-1 interactions with the cells. Firstly, a heterogeneous distribution of IL-1Rs (IL-1R configurations) can be found in several cell types, providing them different capacities to interact with IL-1α and IL-1β (Benjamin and Dower, 1990; Benjamin et al., 1990; Slack et al., 1993). However, how the number of receptors influence the signaling process remains to be elucidated. Secondly, IL-1Rs can be naturally found in their soluble format in several fluids such as blood, peritoneal and synovial fluids (Smith et al., 2003; Michaud et al. 2011; Arend et al., 1994). The soluble receptors interact with both soluble and

2

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

membrane IL-1 ligands and receptors, therefore, inhibiting the signaling process. The physiological role of sIL-1Rs in the control of the inflammation process, as well as their capacity to inhibit IL-1 signaling in different cells, and the relevance of their interactions with the IL-1 ligands and the other receptors have been poorly studied. To shed some light on the open questions regarding the extracellular regulatory pathways of IL-1 signaling and handle the complexity issue, we developed a mathematical model describing the interactions between IL-1 isoforms and the IL-1Rs. A similar model was developed by Kelsey (Kelsey et al., 2008) for the interactions of IL-1β with the IL-1Rs. Unlike Kelsey’s model, we considered the solution and cell membrane as two compartments in which association and dissociation processes are characterized by different kinetic rates. Our model was calibrated with data of IL-1 and IL-1Rs interactions in both solution and cell membrane. With the calibrated model, we simulate the in vitro cells stimulation in different conditions and compute the number of signaling complexes. The model predicts that IL-1α has more ability than IL-1β to form signaling complexes, independently on the cells IL-1Rs configuration. On the other hand, the capacity of the cells to interact with IL-1 mainly depends on the number of IL-1R1 and IL-1R3, and not IL-1R2. The model was extended to study the capacity of soluble forms of IL1Rs (sIL-1R1, sIL-1R2 and sIL-1R3) and IL-1RA to inhibit IL-1 signaling. The extended model predicts that sIL-1R1and IL-1RA have the highest efficiency to block IL-1α ; while sIL-1R2, sIL-1R1and IL1RA have a similar capacity to block IL-1β signaling. Additionally sIL-1R3 has the lowest capacity to prevent IL-1 signaling but it has more specificity to block cells with a low number of IL-1RA. Curiously, the measurements of circulating levels of sIL-1Rs in blood and other fluids suggest that all of them can be as relevant as IL1RA in the control of IL-1 signaling and, therefore the inflammation processes. 2. Methods 2.1. Model equations The model considers that IL-1Rs can move freely across the cell membrane. IL-1 ligands interact with IL-1R1(paths α 1 and β 1) or IL-1R2 (paths α 2 and β 2) from solution, and the resulting complex interacts with IL-1R3 at the cell membrane. This is, by definition, an affinity conversion mechanism (see Fig. 1).

Fig. 1. Scheme for trimeric IL-1-IL-1Rs complexes assembling in the model. Each pathway represents a pair of coupled reactions. The first reaction corresponds to the interaction of soluble ligand with the capturing receptor and the second one to the association of the resulting complex with the IL-1R3. Paths α 1 and α 2 correspond to the assembling of IL-1α -IL-1R1-IL-1R3 and IL-1α -IL-1R2-IL-1R3 complexes, respectively. Paths β 1 and β 2 correspond to the assembling of IL-1β -IL-1R1-IL-1R3 and IL-1β -IL-1R2-IL-1R3 complexes, respectively. Those complexes capable to induce IL-1 signaling are represented in blue. The parameters with unknown value are represented in red. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

The following system of equations describes the reactions represented in paths α 1 and α 2 for the interaction of IL-1α with IL1Rs at the ith cell membrane: i

kα 13 dNα 1 = kα 1 [α ]N1i − k−α 1 Nαi 1 − i Nαi 1 N3i + k−α 13 Nαi 13 dt A

(1)

i

kα 23 dNα 2 = kα 2 [α ]N2i − k−α 2 Nαi 2 − i Nαi 2 N3i + k−α 23 Nαi 23 dt A

(2)

i

dNα 13 kα 13 i i = N N − k−α 13 Nαi 13 dt Ai α 1 3

(3)

i

dNα 23 kα 23 i i = N N − k−α 23 Nαi 23 dt Ai α 2 3

(4)

i i N1i = Ni0 1 − Nα 1 − Nα 13

(5)

i i N2i = Ni0 2 − Nα 2 − Nα 23

(6)

i i N3i = Ni0 3 − Nα 13 − Nα 23

(7)

0 [α ] = [α ] −

 Ncells  i Nα 1 + Niα 2 + Niα 13 + Niα 23 NAV

(8)

An analog system can be written for the interaction with IL-1β by replacing α by β in the subindexes and the concentration of the soluble ligand. Eqs. (1) and (2) describe the dynamics of IL-1α -IL-1R1and IL1α -IL-1R2 complexes respectively. The first and second terms of Eqs. (1) and (2) describe the association and dissociation processes of soluble IL-1α and IL-1R1(IL-1R2). The third and fourth terms of Eqs. (1) and (2) describe the association and dissociation of the complexes IL-1α -IL-1R1(IL-1α -IL-1R2) with the IL-1R3 in the cell membrane. Eqs. (3) and (4) describe the dynamics of IL-1α -IL-1R1-IL-1R3 and IL-1α -IL-1R2-IL-1R3 trimeric complexes due to the association and dissociation of IL-1R3 with IL-1α -IL-1R1 and IL-1α -IL-1R2 complexes respectively. Eqs. (5)–(8) correspond to the conservation law for the mass of IL-1R1, IL-1R2, IL-1R3, and IL-1 molecules. The definition of the model variables and parameters is summarized in Tables 1–3. 2.2. Model parameters The model (see Eqs. (1)–(8)) depends on sixteen kinetic parameters that characterize the interaction between IL-1 ligands and IL1Rs. Eight kinetic rates characterize the interaction of IL-1α and IL-1β with IL-1R1 and IL-1R2 from solution media (kα 1 , k-α 1 , kα 2 , k-α 2 , kβ 1 , k-β 1 , kβ 2 , and k-β 2 ). The value of these parameters was fixed according to reported measurements from surface plasmon resonance (SPR) in previous works in the literature (see Table 2). However, some of the reported values of kinetic rates are inconsistent with the radio-receptor experiments since they correspond to higher affinity than those predicted by the Scatchard plot slopes. Therefore they cannot be used for model fitting (i.e. the parameters measured for the interaction of IL-1β and IL-1R1 and IL-1R2 reported in (Giri et al., 1994)). On the other hand, the kon of IL1α for IL-1R1 (kα 1 ) reported in this paper does not allow a good model fitting due to its low value. As this is the only available report of this parameter, we decide to estimate its value in the model fitting. Eight kinetic rates characterize the interaction of IL-1R3 with the dimeric IL-1α (β )-IL-1R1(2) complexes (kα 13 , k-α 13 , kα 23 , k-α 23 , kβ 13 , k-β 13 , kβ 23 , and k-β 23 ). The value of these parameters cannot be directly measured in the cell membrane scenario, and therefore, is estimated here in model fitting.

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

3

Table 1 Definition of model variables. Symbol

Definition

[α ], [β ]

Concentration of IL-1α and IL-1β , respectively, in the cells containing media

N1i , Ni2 , Ni3

Number of free IL-1R1, IL-1R2, and IL-1R3 on the membrane of ith cell type

Nαi 1 , Niα 2 , Niβ 1 , Niβ 2 Nαi 13 , Niα 23 , Niα 13 , Niα 23

Number of IL-1α -IL-1R1, IL-1α -IL-1R2, IL-1β -IL-1R1, IL-1β -IL-1R2 dimeric complexes on the ith cell type Number of IL-1α -IL-1R1-IL-1R3, IL-1α -IL-1R2-IL-1R3, IL-1β -IL-1R1-IL-1R3, IL-1β -IL-1R2-IL-1R3 trimeric complexes on the ith cell type Table 2 Model parameters whose values were taken from literature. Parameter

Definition

Selected value

k −α 1 kα 2 k −α 2 kβ 1 k −β 1 kβ 2 k −β 2

Dissociation rate of IL-1α from IL-1α -IL-1R1 complex Association rate of IL-1α to IL-1R2 Dissociation rate of IL-1α from IL-1α -IL-1R2 complex Association rate of IL-1β to IL-1R1 Dissociation rate of IL-1β from IL-1β -IL-1R1complex Association rate of IL-1β to IL-1R2 Dissociation rate of IL-1β from IL-1β -IL-1R2 complex

2.2 1.6 1.3 1.5 5.8 1.5 3.1

h−1 × 109 × 102 × 109 h−1 × 109 h−1

(Mh)−1 h−1 (Mh)−1 (Mh)−1

Reference (Giri et al., 1994) (Arend et al., 1994) (Arend et al., 1994) (Smith et al., 2002) (Smith et al., 2002) (Arend et al., 1994) (Arend et al., 1994)

Table 3 Model parameters, whose values are estimated by fitting to data. Parameter

Definition

Range explored in the model fitting

kα 1 kα 13 ∗ k−α 13 ∗ kα 23 ∗ k−α 23 ∗ kβ 13 ∗ k−β 13 ∗ kβ 23 ∗ k−β 23 N1i0

Association rate of soluble IL-1α to IL-1R1 Association rate of free IL-1R3 to IL-1α -ILR1 complex in the cell membrane Dissociation rate of IL-1R3 from IL-1α -ILR1-IL-1R3 complex in the cell membrane Association rate of free IL-1R3 to IL-1α -ILR2 complex in the cell membrane Dissociation rate of IL-1R3 from IL-1α -ILR1-IL-1R3 complex in the cell membrane Association rate of free IL-1R3 to IL-1β -ILR1 complex in the cell membrane Dissociation rate of IL-1R3 from IL-1β -ILR1-IL-1R3 complex in the cell membrane Association rate of free IL-1R3 to IL-1β -ILR2 complex in the cell membrane Dissociation rate of IL-1R3 from IL-1β -ILR1-IL-1R3 complex in the cell membrane Number of IL-1R1 at the membrane of ith cell type in the absence of IL-1

N2i0

Number of IL-1R2 at the membrane of ith cell type in the absence of IL-1

N3i0

Number of IL-1R3 at the membrane of ith cell type in the absence of IL-1

Ai

Relative membrane area with respect to the membrane of B1 cells

(9–180) × 108 (Mh)−1 (1.8 × 10−4 –72) h−1 (3.6 × 10−4 –1.4 × 104 ) (1.8 × 10−4 –72) h−1 (3.6 × 10−4 –1.4 × 104 ) (1.8 × 10−4 –72) h−1 (3.6 × 10−4 –1.4 × 104 ) (1.8 × 10−4 –72) h−1 (3.6 × 10−4 –1.4 × 104 ) A1: 102 –103 A2: 102 –103 A3: 102 –2 × 103 B1: 102 –103 B2: 102 –3 × 103 A1: 103 –5 × 103 A1: 103 –5 × 103 A3: 103 –5 × 103 B1: 0–103 B2: 102 –2 × 102 A1: 102 –8 × 104 A2: 102 –8 × 104 A3: 103 –2 × 104 B1: 10–103 B2: 102 –5 × 103 A1: 0.1–10 A2: 0.1–10 A3: 0.1–10 B2: 0.1–10

∗

h−1 h−1 h−1 h−1

The values of kon rates are normalized by the membrane area of B1 cells.

The system Eqs. (1)–(8) also depends on four cell-specific parameters, which are the number of IL-1Rs and the cell membrane area. The value of these parameters for the five different cells was also estimated in the model fitting. The value of the cell membrane area is always referred to in relative terms to the one of KB cells clone (B1 cells). Two other parameters characterize the experimental conditions: the number of cells and the volume of the cells containing media. The value of these parameters was set to that reported in the original papers (see Table 4). 2.3. Experimental data To properly calibrate the model, it was used published data of both IL-1α and IL-1β binding to human cell lines with different expression of IL-1Rs (Benjamin and Dower, 1990; Benjamin et al., 1990; Slack et al., 1993). In these assays, different concentrations

of IL-1 were added to the cells containing media. After some time, the experiments were stopped, and the number of IL-1 molecules associated to the cells was estimated. The binding assay selected for the model fitting takes four hours after IL-1 addition to the media (except for B1 cells that takes only 2 h). The reported volume and the number of cells in each assay are shown in Table 4.

2.4. Model fitting and parameters estimation In order to reproduce the IL-1 binding assay, the system Eqs. (1)–(8) was numerically solved with the following initial conditions:

Nαi,j1 = Niα,j2 = Niα,j13 = Niα,j23 = 0, [α ] j = [α ]j,0 Nβi,j1 = Niβ,j2 = Niβ,j13 = Niβ,j23 = 0, [β ] = [β ] j

j,0

(9)

4

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113 Table 4 Binding assay characteristics. Cells

Cells origin

Volume

Number of cells

Reference

A1 A2 A3 B1 B2

EBV-transformed cord blood lymphocyte cell lines (clone CB23) EBV-transformed cord blood lymphocyte cell lines (clone CB33) Tumor cell line derived from a patient with African Burkitt’s lymphoma (clone AK778) Human epidermoid carcinoma (clone KB ATCC CCL17) Tumor cell line derived from a patient with AIDS and Burkitt’s lymphoma (PA682BM-2)

300 300 150 150 300

2.7 × 105 2.7 × 105 9 × 106 2.5 × 106 2.7 × 105

(Benjamin and Dower, 1990) (Benjamin and Dower, 1990) (Benjamin et al., 1990) (Slack et al., 1993) (Benjamin and Dower, 1990)

where the superindex j refers to the j-th point in the IL-1 binding assay of the ith cell type. Note that we consider all the cytokine free at the beginning of the experiment. The total amount of IL-1α or IL-1β molecules bound to the ith cell type in the j-th point at time t was computed as:

αˆ bi,j (t ) = Niα,j1 (t ) + Niα,j2 (t ) + Niα,j13 (t ) + Niα,j23 (t ) βˆbi,j (t ) = Niβ,j1 (t ) + Niβ,j2 (t ) + Niβ,j13 (t ) + Niβ,j23 (t )

(10)

To simultaneously fit all the available data, the chi-square function was defined as:



χ2 =

 

Ni N   αˆ bi, j t i − αbi, j

σi,2j

i=1 j=1

i,j

i,j



2 +

Ni N   i=1 j=1

2   βˆbi, j t i − βbi, j σi,2j

(11)

where αb and βb are the experimental value of IL-1α (IL-1β ) bound to ith cell type (N in total) in the j-th point of the assay respectively (Ni in total). Time ti refers to the time consumed in the ith assay. Standard deviations σ α i,j and σ β i,j are considered to be proportional to the mean experimental value (10% due to the good quality of the selected experimental data). An exhaustive search of the chi-square minimums was performed within wide but reasonable ranges for the value of the unknown parameters (see Table 3). The upper limit of association rates in the cell membrane was fixed to the theoretical diffusion limits for the movement of IL-1R chains (see supplementary information 1). The lower limit was taken five orders of magnitude below. The upper limits of the dissociation rates of IL-1R3 and IL1-IL-1R1and IL-1-IL-1R2 complexes were fixed to those estimated in (Wang et al., 2010) for the interaction of sIL-1R3 with IL-1β -IL1R1or IL-1β -IL-1R2 complexes. The lower limit for these rates was fixed in eight orders of magnitude below. The number of IL-1Rs was explored in wider ranges than those estimated in the original works by the Scatchard method. The relative membrane area was explored between 0.1 and 10. The relative membrane area of A1 and A2 cells was assumed the same due to these clones have the same origin and were obtained by the same procedure. The procedure used for model fitting and parameters estimation is similar to that developed in (Ponce et al., 2016). Wolfram Mathematica.11 was used to calculate the chi-square minimums over the parameters values. The common goodness of fit criteria such as Pearson test are not adequate to classify the data coming from the fitting we performed since correlation between parameters are allowed in this process. These correlations make hard to estimate the real degree of freedom of the system. To address this problem, we applied a simple and restrictive criteria based on the distance between the experimental data and the predicted curves. With this criteria, those solutions that predict more than 90 percent of the relative errors within the range (−0.1, 0.1) were chosen as good solutions. Several starting points for minimization method (more than 40 0 0) were randomly generated within the ranges reported in Table 3, obtaining a first set of chi-square minimums. To refine the estimation, a new set of minimums were obtained by generating new starting points (more than 10 0 0) around the best solution

μL μL μL μL μL

within the first set of minimums (within +/− 20% of the value of each parameter). A deeper exploration of the confidence interval of the fitting was performed for the values of the kinetic coefficients in the cell membrane. For that, a new set of solutions were generated by varying the value of these parameters by pairs of kon and koff (104 iterations by pair). In this procedure, the values of the remaining parameters were fixed to those obtained in the original selected solution. This procedure generated a new set of solutions that do not correspond to chi-square local minimums and was filtrated by the criteria explained above. To test how the variations in the parameters whose values were taken from the literature, it was performed a robustness analysis (Supplementary Information 4). In this analysis, the values of the kinetic rates for the interaction of IL-1 with IL-1R1 and IL-1R2 were varied in a range between 0.5 and 1.2 times the value used in the model fitting. The comparison of both the cell-dependent parameters and the membrane kinetic coefficients shows that the variations in the kinetic rates of solution interaction can be compensated with small variations in the other parameters. In this way, both the estimated range for the parameters values and the relation between them are conserved. 3. Results and discussion 3.1. The model properly fits the experimental data To study the extracellular regulation of IL1 we formulated a mathematical model that describes the kinetics of interaction of IL-1α or IL-1β with cells expressing IL-1Rs. The model was calibrated with data of IL-1α and IL-1β binding to cells with different IL-1Rs configurations. As expected, the model can fit each data set individually in a wide range of parameter values (results not shown). Interestingly, a successful simultaneous model fitting was obtained using the five available data sets. An example of a good simultaneous model fitting is shown in Fig. 2. The histogram of relative error frequencies (last panel) shows a distribution similar to a normal distribution, where more than 90 percent of the values lay within −0.1 and 0.1, confirming the quality of the fitting. 3.2. Estimated values of cell-dependent parameters Fig. 3 shows the estimated value of the number of IL-1Rs and the membrane area for all cells used in the model fitting and the respective chi-square value. Each point in the graph corresponds to a good solution selected by the criteria explained in Section 2. A low variation of the estimated value between different solutions was obtained for most of these parameters (i.e. the number of IL1R1and IL-1R2 in A1 cells, see Table 5), meaning that they can be considered as structural identifiable according to Raue’s definition (Raue et al., 2009). This result can be explained by the fact that the experimental data contains restraining information of two different molecules interacting with the same IL-1Rs but with different affinities. The sensitivity analysis (Supplementary Information 4) confirms that varying the value of the parameters in more than

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

5

Fig. 2. The model simultaneously fits the experimental data. Dots correspond to the experimental data from (Benjamin and Dower, 1990; Benjamin et al., 1990; Slack et al., 1993) and solid lines to the predicted Scatchard plot with the fitted model. Labels above each curve indicate the corresponding ligand IL-1α or IL-1β . Label on the top right of each graph indicates the cell type. The parameter values correspond to the best solution (last four columns of Table 5 and last column of Table 6). Last graphic corresponds to the histogram of relative errors between model prediction and experimental data.

20% of its value makes the solutions incompatible with the experimental data. In the case of cells with excess of IL-1R3 (A1 and A2 cells), both the number of IL-1R3 and the relative membrane area are practically non-identifiable because of how these parameters appear in the model, in spite of the low variation obtained in the fitting process.

The estimated configurations of IL-1Rs predict two groups of cells: one with excess of IL-1R2 over IL-1R1(cells A1-A3); and other with the inverse proportion of these IL-1Rs (B1 and B2 cells). Within the group A, it was estimated that IL-1R3 is in excess with respect to IL-1R1 and IL-1R2. On the other hand, the estimated number of IL-1R3 was slightly lower than IL-1R1 in B1 cells and

Table 5 Estimation of cell-dependent parameters. Estimated range

Best Fit

cell

ILR1

IL-1R2

IL-1R3

A

IL-1R1

IL-1R2

IL-1R3

A

A1 A2 A3 B1 B2

136–154 151–177 358–428 1060–1192 359–374

2292–2482 2042–2361 16878–17252 – 18–35

6514–30331 9674–46193 9445–54009 911–2763 204–258

0.08–0.26 0.08–0.26 0.3–0.9 1 1.1–6.0

140 155 382 1066 365

2388 2220 17046 0 29

17453 27580 31319 2573 229

0.25 0.25 0.87 1 5.1

6

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

Fig. 3. Estimation of the number of IL-1Rs and membrane area in different cells. Panels A1-B2 show the value of the number of IL-1Rs in each local minimum of chi-square function selected with the criteria explained in the methods section. Each point represents the chi-square (y axis) and the value of the parameter for one solution (x axis). In panels A1-B2, the blue, yellow and green colors correspond to the number of IL-1R1, IL-1R2 and IL-1R3 respectively. Last panel shows the estimated value of the relative membrane area with respect to B1 cells. Each color represents the estimation for a specific cell, which is labeled in the top of the graph. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

similar in B2 cells. Due to the size of the data sample, it is not possible to get any important conclusion about the preferential IL1Rs configuration depending on the cell type. The estimation of the relative membrane area (last panel in Fig. 3 and last column of Table 3) shows the lowest values in A1 and A2 cells. This result was expected since these cells were obtained from cord blood lymphocytes, which are described as small size cells (less than 6 μm of diameter (Rodrigues et al., 2010). However, different areas were estimated for cells with similar origin (example A3 and B2 cells derived from Burkitt’s lymphomas). We interpreted that the predicted area in these cases could be related to the effective area for the movement of IL-1Rs more than the size of the cell membrane.

3.3. Estimation of kinetic coefficients at the cell membrane Fig. 4 shows the estimated values for kinetic coefficients characterizing the interactions at the cell membrane. Each point corresponds to a pair of kon and koff within a selected solution. The estimated values of these rates vary in more than 3 orders of magnitude (see Table 6), indicating that these are practical nonidentifiable parameters (Raue et al., 2009). However, an accurate estimation was obtained for the ratio kon/koff in the reactions between IL-1R3 and IL-1α -IL-1R1, IL-1α -IL-1R2 and IL-1β -IL-1R2 complexes, which is by definition, the affinity of these reactions. For the interaction between IL-1R3 and IL-1β -IL-1R2 complex, an inaccurate estimation was obtained of both the affinity and ki-

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

7

Fig. 4. Estimation of kinetic coefficients for the interaction of IL-1-IL-1R1(2) complexes with IL-1R3 at the cell membrane. In the top of graphs is labeled the corresponding ligand. Blue and yellow colors correspond to the interaction of IL-1R3 with the complex IL-1-IL-1R1and IL-1-IL-1R2 respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

netic coefficients (see Fig. 4, panel b in yellow, and Table 6). In this case, only a maximum value for the affinity was obtained. This maximum value is more than 2 orders of magnitude lower than the affinity predicted for the interaction of IL-1R3 with IL-1β -IL1R. Interestingly, this difference was also measured for the interaction of soluble IL-1R3 and the complexes IL-1β -IL-1R1and IL-1β IL-1R2 by SPR (Wang et al., 2010). The previous results indicate that affinities, rather than kinetic coefficients, seem to be constrained in the model fitting. This result is explained by the fact that the radio-receptor assays are typically experiments on equilibrium conditions. Therefore, it must be expected that experimental data is less informative regarding the kinetic to reach such equilibrium. The kon rate of interaction between IL-1α with IL-1R1 was precisely estimated between 2.1 to 2.4 times the one reported by Giri et al. (1994). 3.4. Model predictions: cells are more sensitive to IL-1α than IL-1β The calibrated model was used to study how IL-1-IL-1Rs system output is regulated by ligand concentrations and the configuration of IL-1Rs. First, we analyzed the possible differences of IL-1α or IL1β in the formation of IL-1 signaling complexes. This analysis was performed for a high number of possible configurations of IL-1Rs, including those estimated for A1-B2 cells. As it is shown in Fig. 5, IL-1α has a greater capacity to form signaling complexes than IL-1β . In all analyzed configurations, the Ec50 value is lower for IL-1α , which is mainly due to the higher affinity of IL-1α for IL-1R1. Additionally, for cells in which the

Table 6 Estimation of kinetic coefficients. Parameter

Estimated range

Best fit

kα 1 kα 13 k-α 13 Kα 13 = kα 23 k-α 23 Kα 23 = kβ 13 k-β 13 Kβ 13 = kβ 23 k-β 23 Kβ 23 =

(1.9–2.2) × 109 (Mh) −1 (7.5 × 10−2 - 70)h−1 (5.8 × 10−1 −1.1 × 103 )h−1 3.7 × 10−2 −1.7 × 10−1 (3.5 × 10−4 - 36)h−1 (3.8 × 10−3 −1.0 × 104 )h−1 1.2 × 10−3 −2.8 × 10−1 (5.9 × 10−3 −69)h−1 (9.6 × 10−2 −2.9 × 103 )h−1 1.6 × 10−2 −6.1 × 10−2 (1.8 × 10−4 −1.7)h−1 (9.6 × 10−2 −1.0 × 104 )h−1 1.4 × 10−2 −2.9 × 10−4

2.0 1.4 1.1 1.3 6.5 5.0 1.3 9.7 2.6 3.6 5.8 1.3 4.3

kα 13 / k-α 13

kα 23 / k-α 23

kβ 13 / k-β 13

kβ 23 / k-β 23

× 109 (Mh)−1 × 10 h−1 × 102 h−1 × 10−1 × 10−2 h−1 h−1 × 10−2 h−1 × 102 h−1 × 10−2 × 10−1 h−1 × 103 h−1 × 10−4

Fig. 5. IL-1α is more active than IL-1β . Graphs show the predicted number of signaling complexes in equilibrium as a function of the concentration of IL-1α (blue) or IL-1β (yellow). In the top of the graphs is labeled the cells types from which cellspecific parameters were selected. The model parameters correspond to the best fit (last 4 columns of Table 5 and last column of Table 6). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

number of IL-1R1or IL-1R3 are not in clear excess to each other, it was predicted that the signaling plateau is higher for IL-1α than of IL-1β (see Fig. 5, B1 cells). This result can be explained by the fact that the complex IL-1α -IL-1R1 has a higher affinity for IL-1R3 in the cell membrane than the complex IL-1β -IL-1R1. The differences between the response of cells to IL-1α and IL1β stimulation have been poorly studied. However, in (Neumann et al., 20 0 0), it was obtained that IL-1α has more capacity than IL-1β to induce IL6 secretion in HaCaT cells, which is in agreement with our prediction. In another system (different to the human one), in the work of Tanikawa et al. (2008), it was found that IL-1α has a larger capacity than IL-1β to stimulate prostaglandin synthesis in endometrial stromal cells. The previous result alerts about the role of IL-1α in physiological and inflammatory conditions. Most of the published works regarding IL-1 only refers to IL-1β and not IL-1α , probably due to the higher levels of IL-1β observed in some scenarios. However, our prediction indicates that IL-1α can be as relevant as IL-1β , even at lower concentrations. For example, in (Aguilar-Santelises et al., 1992) it was measured around the double of IL-1β serum concentration than IL-1α in healthy donors (14 pg/mL and 26 pg/mL, respectively), and similar concentration of both ligands in patients with monoclonal lymphocytosis (18 pg/mL and 23 pg/mL, respectively) or progressive B-lymphocytic leukemia (19 pg/mL and 22 pg/mL, respectively). Those serum concentrations are low compared with IL-1α and IL-1β Ec50 values, but higher concentrations are expected at the tissue level and inflammation sites. Considering that serum measurements reflect somehow the dynamic of the tissues, it is expected that the concentration of IL-1α will be around

8

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

the half of IL-1β in the normal or inflamed tissues. According to the model, in these scenarios, IL-1α signaling could be similarly, or even more relevant than the one of IL-1β , for cells with any IL-1Rs configurations. Until now, there is no clear evidence that the difference between IL-1α and IL-1β interactions with membrane IL-1Rs impacts the quality of IL-1 signaling. However, recent works suggest that IL-1α and IL-1β are involved in different immunological processes (Lee et al., 2011; Rider et al., 2011). In the Rider’s work, it was found that IL-1α stimulates different myeloid cells than IL-1β . They notice that IL-1α is accumulated in the early phase of inflammation and recruits mostly neutrophils while IL-1β is more related to the later recruitment of macrophages. The results of Lee reinforce the idea of the involvement of IL-1α in the recruitment of neutrophils, in this case in the chronic inflammation induced by hydrocarbon oil.

3.5. Model predictions: modulation of IL-1 signaling by IL-1Rs expression The quantitative contribution of IL-1Rs to cell signaling capacity and sensibility have been poorly studied experimentally. Here, we used the calibrated model to predict how the IL-1Rs configuration influences the capacity of cells to interact with IL-1α and IL-1β . For that, we computed the dependence of signaling plateau and EC50 value with the number of a given IL-1R type, while the others IL-1Rs and the relative membrane area remain constant. This procedure was performed for the three types of IL-1Rs (see Fig. 6 for A1 and B1 cells). As expected, the model predicts a relevant role of IL-1R1and IL1R3 in controlling IL-1 signaling. However, the contribution of both receptors varies depending on IL-1Rs configurations. For example, for A1 cells it was predicted that EC50 value decreases with the number of IL-1R3 and increases slightly with the number of IL1R1, for both ligands (see panels A and C in Fig. 6). The signaling plateau in this cell type increases mostly with the number of IL1R1 (panels E and G of Fig. 6). On the other hand, for B1 cells, the EC50 value decreases with the number of both IL-1R1 and IL1R3 but shows a higher dependency on the number of IL-1R1(see panels B and D in Fig. 6). Moreover, the signaling plateau increases with both IL-1R1and IL-1R3 and reaches a maximum value, showing a larger variation with the number of IL-1R3 (see panels F and H in Fig. 6). Interestingly, EC50 value and the signaling plateau are independent or slightly decrease with IL-1R2. This behavior is due to the low affinity of IL-1-IL-1R2 complex with the IL-1R3. Therefore, it is expected to be more relevant in cells with a low number of IL-1R3 (like B1 cells). Overall, these results indicate that, at the cell level, IL-1 signaling is mainly controlled by the expression of IL-1R1and IL-1R3, and not IL-1R2. This result is consistent with the predictions of Kelsey’s model (Kelsey et al. 2008). The low or almost null capacity of the IL-1R2 to prevent IL-1 signaling was experimentally tested by Neumann et al. (20 0 0). In this work, it was found that the capacity of IL-1α and IL-1β to stimulate HaCaT cells was the same when the cells are forced to express a large number of IL-1R2 or a mutated IL-1R2 isoform that is unable to form IL-1-IL-1R2-IL-1R3 complexes. This result reinforces the idea that the interaction of IL-1-IL-1R2 with the IL-1R3 is not important for the controlling of IL-1 signaling. As the number of IL-1R3 determines the capacity of the cells to interact with IL-1, any molecule that share this correceptor can regulate this interaction. Therefore, the capacity of the cells to interact with IL-1 in vivo will be determined indirectly by the expression of ST2 and IL-36Rα and the presence in the media of IL33 and any isoform of IL-36.

Fig. 6. IL-1Rs expression modulates the capacity of cells to interact with IL-1. Panels A-D show the dependency of EC50 value with the number of IL-1R1(blue), IL1R2 (yellow) and IL-1R3 (green). Panels E-H show the dependency of the maximum number of signaling complexes with the number of IL-1Rs (the same color code). Panels A, B, E and F correspond to predictions for IL-1α ; and panels C, D, G and H show the predictions for IL-1β . The configurations estimated for A1 and B1 cells were used to obtain the graphics in the left and right columns respectivelly. Filled bigger circles correspond to the original numbers of IL-1Rs. The parameters used in the simulations correspond to the best solution of model fitting (the last four columns of Table 5 and last column of Table 6). The end time of these simulations were 10 h to unsure the equilibrium of the system. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

3.6. Model predictions: modulation of IL-1 signaling by soluble form of IL-1Rs and IL-1RA IL-1Rs can be found in soluble form in normal conditions. In principle, the soluble form IL-1Rs in the extracellular media prevents the formation of signaling complexes by different mechanisms and with different efficiency. To explore how sIL-1Rs regulate IL-1α and IL-1β interaction with the cells, the model (Eqs. (1)– (8)) was extended to simulate the scenario in which the solution media has a certain concentration of IL-1α or IL-1β , sIL-1Rs, and IL-1RA (see the extended model in Supplementary Information 2 and 3). Fig. 7 shows the number of signaling complexes formed with IL-1α or IL-1β in two different cells, for increasing concentrations of sIL-1Rs and IL-1RA. In these simulations, it was considered a fixed initial concentration of free IL-1α or IL-1β corresponding to the EC50 value. We obtain that sIL-1R1 has the greatest capacity to prevent IL-1α signaling for both cell configurations followed by IL-1RA. On the other hand, sIL-1R1, sIL-1R2 and IL-1RA have a similar capacity to prevent IL-1β signaling. Moreover, sIL-1R3 has the lowest capacity to prevent IL-1 signaling. However, sIL-1R3 is the only sIL-1R that preferentially inhibits certain type of cells: those cells with a lower number of IL-1R3.

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

Fig. 7. Antagonist capacity of soluble form of IL-1Rs and IL-1RA. Graphs show the number of signaling complexes when the cells are stimulated with IL-1α (panels A and B) or IL-1β (C and D) at a concentration corresponding to EC50 value. Increasing concentration of sIL-1Rs (sIL-1R1in blue, sIL-1R2 in yellow and sIL-1R3 in green) and IL-1RA (red) are added to the media. In panels A and C it was used the IL-1Rs configuration estimated for A1 cells, and, in panels B and D, the configuration estimated for B1 cells. The model parameters correspond to the best solution (see last four columns of Table 5 and last column of Table 6). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

The capacity of IL-1RA, sIL-1R1 and sIL-1R2 to inhibit IL-1 signaling shows an unexpected independence with the cell configuration. Understandably, IL-1RA does not have any cell specificity since it only prevents IL-1 signaling by interacting with IL-1R1 at the cell membrane. The same behavior will show all the molecules that prevent IL1 signaling by using this mechanism. For example, IL-38 also binds to IL-1R1 although with more than 10 times lower affinity than IL-1RA (Lin et al., 2001). Until now, there is no evidence of the recruitment of IL-1R3 by IL-38, thus, it should act in vivo as an IL-1 antagonist and display similar independence on the cell configuration to prevent IL-1 signaling. The independency of the capacity of sIL-1R1and sIL-1R2 to prevent IL1 signaling on the cell receptors configuration is not obvious since they can interact with IL-1R3 once they capture IL-1 in solution, affecting the availability of IL-1R3 in the cell membrane. To evaluate the contribution of this interaction, we performed similar simulations but setting its affinity to zero. The output of these simulations was very similar to those of the previous one (result not shown), indicating that soluble IL-1-sIL-1Rs complexes do not compete with membrane IL-1-IL-1R1for IL-1R3, even in the cells with a low number of IL-1R3. Therefore, the main mechanism of IL-1 signaling control by sIL-1R1 and IL-1R2 is the IL-1 sequestration in solution, which is independent on the cell configuration. The ability of membrane IL-1R3 to interact with IL1α (β )-IL-1R1 complex also explains the low capacity of sIL-1R3 to prevent IL-1 signaling. The capacity of soluble receptors has been tested in vitro in several works (Netea et al., 1999; Kollewe et al. 20 0 0; Smeets et al., 2005). In the work of Kollewe et al. (20 0 0) it was found that sIL1R2 is highly efficient to blockade IL-1β but not IL-1α , which is consistent with its affinity for these molecules. On the other hand, in the work of Smeets et al. (2005) it was found that sIL-1R3 is more capable to prevent the IL-1 signaling in cells with a low number of IL-1R3, which is in agreement with our predictions. Moreover, the concentrations of sIL-1R3 to achieve the significant blockade of IL-1 was found to be higher in comparison to those used for sIL-1R1and sIL-1R2 (μg/mL vs ng/mL), which is also consistent with our predictions.

9

The biological relevance of the soluble IL-1Rs depends on their physiological concentration. The serum concentrations of sIL-1Rs and IL-1RA (Meier et al., 20 0 0; Smith et al., 2003; Alshevskaya et al., 2015) are low compared to the predicted IC50 values (necessary concentration to get 50% of inhibition). Interestingly, the serum concentrations correlates with the respective IC50 value, which suggests that all sIL-1Rs are similarly important for the IL1 signaling regulation. Given the physiological relevance of IL-1RA, it can be inferred that sIL-1Rs also play a significant role in the control of inflammatory processes. The results of the previous section indicate that IL-1R2 does not control IL-1β signaling at the cell level. On the other hand, in this section we predicted that IL-1R2 is capable to prevent IL-1β signaling with high efficiency in the soluble format. Those results suggest that the role of IL-1R2 is the control of the availability of IL1β either in soluble form or a paracrine manner around the cells expressing it. The IL-1β sequestration by IL-1R2 can be a mechanism of action of cells with regulatory phenotype such as Regulatory T cells (Mercer et al., 2010) and M2 polarized macrophages (Mia et al., 2014); and a mechanism of immune suppression of glucocorticoids (Varga et al., 2016). 4. Concluding remarks Due to the capacity of IL-1 to activate the immune system, several mechanisms of control have been established during the evolution. These mechanisms have both intracellular and extracellular scope. The intracellular mechanisms for controlling IL-1 signaling cascade and secretion have been well studied experimentally in spite of their complexity (Weber et al., 2010). On the other hand, several extracellular mechanisms such as decoys, soluble receptors and antagonist ligands have been described for IL-1 system that control the formation of IL-1 signaling complexes in the cells membrane. However, little is know about their relevance and limits in different scenarios. We developed a mathematical model to study the possible extracellular regulation of IL-1 signaling by the components of IL-1 system. The model was calibrated with experimental data of IL1α and IL-1β binding to cells with different configurations of IL1Rs. The model fitting allowed us to estimate, for the first time, the affinity of IL-1-IL-1R1and IL-1-IL-1R2 complexes with IL-1R3 in the cell membrane. These parameters are difficult to determine by the common techniques used to study protein-protein interactions. Moreover, these parameters are useful to simulate longer time processes in which IL-1 is involved such as tumor or tissue growing and vascularization, immune response development, cell migration. Additionally, several aspects of IL-1 signaling regulation were deduced from the model that are expected to be valid for these scenarios: •

•

•

•

•

Independently on the IL-1Rs configuration, cells are more sensitive (lower Ec50 value) to IL-1α than IL-1β stimulation. Moreover, cells with a comparable number of IL-1R1and IL-1R3 have higher responsiveness capacity (signaling plateau) to IL-1α than IL-1β . IL-1 signaling at the cell level is mainly controlled by the expression of IL-1R1and IL-1R3. The role of IL-1R2 is the control of the availability of IL-1β either as soluble format or in a paracrine way around the cells expressing it. sIL-1R1 is the best controller of IL-1α signaling followed by IL1RA; while sIL-1R2, sIL-1R1 and IL-1RA have a similar capacity to control IL-1β . sIL-1R3 has the lowest capacity to control IL-1 signaling among the sIL-1Rs, but is the only one that has specificity for those cells with a low number of IL-1R3.

10

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113

Several experimental works have somehow tested our predictions. However, a deeper exploration would be necessary to generalize these results. On one side, the capacity of IL-1α and IL-1β to stimulate cells and the capacity of soluble IL-1Rs to prevent IL-1 signaling can be tested by measuring any component of the signaling cascade such as the NfkB activation; or an effect such as IL6 secretion. We consider that it would be also interesting to test whether the signaling cascade of IL-1α differs qualitatively from the one induced by IL-1β . In the works of Lee and Rider (Lee et al., 2011; Rider et al., 2011) it was found that IL-1α is more involved in the recruitment of neutrophils while IL-1β is more related to the development of macrophages. In the work of Trebec-Reynolds (Trebec-Reynolds et al., 2010) it was found that large osteoclasts respond differently to IL-1α and IL-1β stimulation. Our results indicate that the complex IL-1α -IL-1R1has a higher affinity for the IL-1R3 than the complex IL-1β -IL-1R1. Recent results pointed out that the stability of the aggregating complexes at the cell membrane determines the signaling kinetics (Freed et al., 2017). According to Freed’s results, the affinity for IL-1R3 might be associated with differences in the signaling cascade of IL-1α and IL-1β . The influence of the number of IL-1Rs (specially IL-1R2) on the capacity of the cell to interact with IL-1 (Ec50 value) can also be tested in vitro especially in cells population with different states of differentiation. The cell-to-cell variability analysis performed by Cotari et al. (2013) for the IL-2 system would be an interesting approach to this problem. It is also remarkable from the literature that with the knowledge of the antagonist capacity of soluble IL-1 receptors, little is written about their role in physiological conditions and inflammatory diseases. It is also unclear what is the role of the membrane attached IL-1R2 and how different is from the soluble format. On the contrary, the IL-1RA has been associated with the control of inflammation in several tissues (Arend and Gabay, 20 0 0). Only some works suggest that sIL-1R2 helps to control inflammatory diseases such as arthritis since the level of sIL-1R2 negatively correlates with the disease score (Peters et al., 2013). According to our predictions about the capacity of soluble IL-1Rs to prevent IL-1 signaling and the measurement of their blood levels, it would be expected that they are as important as IL-1RA in regulating the IL-1 signaling. Following this line of thought, it could be inferred that sIL-1R3 is involved in the control of IL-1 related inflammation in the endometrium since its level decreases with the development of endometriosis (Michaud et al., 2011). At the same time, an imbalance in the level of the soluble receptors could also be a cause of the development of inflammatory diseases. In the modeling of the interaction between IL-1 and the cells it was necessary to make some assumptions and simplifications that can affect model predictions. For example, we considered, for simplicity, a homogeneous distribution. It is possible that a nonhomogeneous spatial distribution of IL-1Rs across the membrane significantly modifies the way in which they interact. On the other hand, IL-1Rs are expressed in many cell types and therefore, kinetic coefficients may vary depending on the membrane composition. Additionally the kinetic rates used in the model for the interaction of IL-1 and IL-1Rs in solution were assumed constants. However, they were measured in a very simple composition media. A more complex media (i.e. the blood, or synovial fluid) with different compositions (different pH values as well) may also influence the way in which molecules interact with each other by affecting the kinetic rates. Processes such as internalization and receptors expression are not considered in this work. However, for similar systems, we found that these processes do not significantly affect the model calibration-prediction, especially for large time compared to one needed to reach the equilibrium point. Other possible interac-

tions that are not considered in the model may also affect our predictions. For example, in the structure of IL-1β -IL-1R1-IL-1R3 and IL-1β -IL-1R2-IL-1R3 complexes it was found a surface of interaction between IL-1R1 and IL-1R2 with IL-1R3 (Wang et al., 2010; Thomas et al., 2012) which suggest that it is possible certain pre-association between these IL-1Rs. This possible pre-association processes may also affect both the model predictions and the interpretation of the experimental results. Authors contribution LFP, KL and KG-M developed the mathematical model and generate the idea of the work. LFP performed the data searching and processing and the numerical solution of equations. LFP KG-M and KL wrote the paper and analyzed the results. Declaration of Competing Interest All authors concur with the submission. All the persons participating in the work are included in the list of authors. Non-funding were used for the studies in the manuscript. The manuscript does not contain experiments using animals. The manuscript does not contain human studies. There isn’t any financial/commercial conflicts of interests. Acknowledgments We thanks to Dra Marina Ferreira from Imperial College of London, Ms Alina and Robert Lawlor from Becton & Dickinson in Ireland, and MSc Mary Kudarybergenova from the Center for Molecular Simulations at the University of Calgary for the writing assistance and proof reading of the manuscript configuration. Supplementary materials Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.jtbi.2019.110113. References Aguilar-Santelises, M., et al., 1992. Serum levels of helper factors (IL-1α IL-1β and IL-6), T-cell products (sCD4 and sCD8), sIL-2R and b2-microglobulin in patients with B-CLL and benign B lymphocytosis. Leuk. Res. 16 (6–7), 607–613. doi:10. 1016/0145-2126(92)90 0 09-V. Alshevskaya, A.A., et al., 2015. Differences of IL-1β receptors expression by immunocompetent cells subsets in rheumatoid arthritis. Mediators Inflamm. 2015. doi:10.1155/2015/948393. Arend, W.P., Gabay, C., 20 0 0. Physiologic role of interleukin-1 receptor antagonist. Arthritis Res. 2 (4), 245–248. doi:10.1186/ar94. Arend, W.P., et al., 1994. ‘Binding of IL-1ɑ, IL-1β , and IL-1 receptor antagonist by soluble IL-1 receptors and levels of soluble IL-1 receptors in synovial fluids. J. Immunol. 153 (10), 4766–4774 PMID: 7963543. Banerjee, M., Saxena, M., 2012. Interleukin-1 (IL-1) family of cytokines: role in type 2 diabetes. Clinica Chimica Acta 1163–1170. doi:10.1016/j.cca.2012.03.021. Benjamin, D., Dower, S.K., 1990. Human B cells express two types of interleukin-1 receptors. Blood 75 (10), 2017–2023. doi:10.1182/blood.V75.10.2017.2017. Benjamin, D., et al., 1990. Heterogeneity in interleukin (IL)-1 receptors expressed on human b cell lines. J. Biol. Chem. 265 (17), 9943–9951. Chen, M.-F., et al., 2012. Role of interleukin 1 beta in esophageal squamous cell carcinoma. J. Mol. Med. 90 (1), 89–100. doi:10.10 07/s0 0109-011-0809-4. Cotari, J.W., et al., 2013. Cell-to-Cell variability analysis dissects the plasticity of signaling of common ɣ chain cytokines in t cells. Sci. Signal 6 (266), 1–10. doi:10.1126/scisignal.2003240. Cullinan, E.B., et al., 1998. IL-1 receptor accessory protein is an essential component of the IL-1 receptor. J. Immunol. (Baltimore, Md. : 1950) 161 (10), 5614–5620. Dinarello, C.A., 2009. Immunological and inflammatory functions of the interleukin1 family. Annu. Rev. Immunol. 27 (1), 519–550. doi:10.1146/annurev.immunol. 021908.132612. Dinarello, C.A., 2010. Why not treat human cancer with interleukin-1 blockade? Cancer Metastasis Rev. 317–329. doi:10.1007/s10555-010-9229-0. Dinarello, C.A., van der Meer, J.W.M., 2013. Treating inflammation by blocking interleukin-1 in humans. Semin. Immunol. 469–484. doi:10.1016/j.smim.2013.10. 008. Freed, D.M., et al., 2017. ‘EGFR ligands differentially stabilize receptor dimers to specify signaling kinetics. Cell 171, 1–13. doi:10.1016/j.cell.2017.09.017.

L.F. Ponce, K. García-Martínez and K. León / Journal of Theoretical Biology 487 (2019) 110113 Garlanda, C., et al., 2013. Decoys and regulatory “receptors” of the il-1/toll-like receptor superfamily. Front. Immunol. 4 (180), 1–13. doi:10.3389/fimmu.2013. 00180. Giri, J.G., et al., 1994. Elevated levels of shed type II IL-1 receptor in sepsis. Potential role for type II receptor in regulation of IL-1 responses. J. Immunol. 153 (12), 5802–5809. Hagemann, T., et al., 2008. Re-educating tumor-associated macrophages by targeting NF-κ B’. J. Exp. Med. 205 (6), 1261–1268. doi:10.1084/jem.20080108. Kay, J., Calaberase, L., 2004. The role of interleukin-1 in the pathogenesis of rheumatoid arthritis. Rheumatology 43, 1–8. doi:10.1093/rheumatology/keh201. Kelsey, J., et al., 2008. ‘A stochastic model of the interleukin (IL)-1β network. Artific. Immune Syst. Springer 1–11. doi:10.1007/11536444. Kollewe, C., et al., 20 0 0. The first two N-terminal immunoglobulin-like domains of soluble human IL-1 receptor type II are sufficient to bind and neutralize IL-1β . FEBS Lett. 487, 189–193. doi:10.1016/s0014- 5793(00)02345- 0. Lee, P.Y., et al., 2011. IL-1 α modulates neutrophil recruitment in chronic inflammation induced by hydrocarbon oil. J. Immunol. 186 (3), 747–1754. doi:10.4049/ jimmunol.1001328. Lin, H., et al., 2001. Cloning and characterization of IL-1HY2, a novel interleukin1 family member. J. Biol. Chem. 276 (23), 20597–20602. doi:10.1074/jbc. M010 09520 0. Meier, C.A., et al., 20 0 0. IL-1 receptor antagonist serum levels are increased in human obesity: a possible link to the resistance to leptin? J. Clinic. Endocrinol. Metabol. 87 (3), 1184–1188. doi:10.1210/jcem.87.3.8351. Mercer, F., et al., 2010. Expression and function of TNF and IL-1 receptors on human regulatory T cells. PLoS ONE 5 (1), 1–12. doi:10.1371/journal.pone.0 0 08639. Mia, S., et al., 2014. An optimized protocol for human M2 macrophages using M-CSF and IL-4/IL-10/TGF- b Yields a dominant immunosuppressive phenotype. Scand. J. Immunol. 79 (5), 305–314. doi:10.1111/sji.12162. Michaud, N., et al., 2011. Decreased concentrations of soluble interleukin-1 receptor accessory protein levels in the peritoneal fluid of women with endometriosis. J. Reprod. Immunol. 92, 69–73. doi:10.1016/j.jri.2011.08.001. Netea, G.M., et al., 1999. ‘Soluble murine IL-1 receptor type I induces release of constitutive IL-1α . The J. Immunol. 162, 4876–4881. Neumann, D., et al., 20 0 0. The membrane form of the type II IL-1 receptor accounts for inhibitory function. J. Immunol. 165 (6), 3350–3357. doi:10.4049/jimmunol. 165.6.3350. Peters, V.A., et al., 2013. IL-1 receptor 2 (IL-1R2) and its role in immune regulation. Brain Behav. Immun. 32, 1–8. doi:10.1016/j.bbi.2012.11.006. Ponce, L.F., et al., 2016. Quantitative contribution of IL2Rγ to the dynamic formation of IL2-IL2R complexes. PLoS ONE 11 (5). doi:10.1371/journal.pone.0155684. Pope, R.M., Tschopp, J., 2007. The role of interleukin-1 and the inflammasome in gout: implications for therapy. Arthritis Rheum. 56 (10), 3183–3188. doi:10. 1002/art.22938. Raue, A, et al., 2009. Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood. Bioinformatics 25 (15), 1923–1929. doi:10.1093/bioinformatics/btp358. Ricote, M., et al., 2004. Interleukin-1 (IL-1α and IL-1β ) and its receptors (IL-1R, IL-1R2, and IL-1RA) in prostate carcinoma. Cancer 100 (7), 1388–1396. doi:10. 1002/cncr.20142.

11

Rider, P., et al., 2011. IL-1 α and IL-1 β recruit different myeloid cells and promote different stages of sterile inflammation. J. Immunol. 187 (9), 4835–4843. doi:10. 4049/jimmunol.1102048. Rodrigues, G.C, et al., 2010. Ultrastructural characterization of bovine umbilical cord blood cells 1. Pesquisa Veterinária Brasileira 30 (10), 897–902. doi:10.1590/ S0100-736X2010 0 010 0 0 015. Saijo, Y., et al., 2002. Proinflammatory cytokine IL-1 promotes tumor growth of lewis lung carcinoma by induction of angiogenic factors: in vivo analysis of tumor-stromal interaction. J. Immunol. 169 (1), 469–475. doi:10.4049/jimmunol. 169.1.469. Shaftel, S.S., et al., 2008. The role of interleukin-1 in neuroinflammation and Alzheimer disease: an evolving perspective. J. Neuroinflammat. 5 (1), 1–12. doi:10.1186/1742- 2094- 5- 7. Slack, J., et al., 1993. Independent binding of interleukin-1 alpha and interleukin-1 beta to type I and type II interleukin-1 receptors. J. Biol. Chem. 268 (4), 2513–2524. Smeets, R.L., et al., 2005. Soluble interleukin-1 receptor accessory protein ameliorates collagen-induced arthritis by a different mode of action from that of interleukin-1 receptor antagonist. Arthritis Rheum. 52 (7), 2202–2211. doi:10. 1002/art.21108. Smith, D.E., et al., 2002. A single amino acid difference between human and monkey interleukin (IL)-1β dictates effective binding to soluble type II IL-1 receptor. J. Biol. Chem. 277 (49), 47619–47625. doi:10.1074/jbc.M206636200. Smith, D.E., et al., 2003. The soluble form of IL-1 receptor accessory protein enhances the ability of soluble type ii IL-1 receptor to inhibit IL-1 action. Immunity 18, 87–96. doi:10.1016/s1074- 7613(02)00514- 9. Tanikawa, M., et al., 2008. ‘Regulation of prostaglandin biosynthesis by interleukin-1 in cultured bovine endometrial cells. J. Endocrinol. 199, 425–434. doi:10.1677/ JOE- 08- 0237. Thomas, C., et al., 2012. Structure of the activating IL-1 receptor signaling complex. Nat. Struct. Mol. Biol. 19 (4), 455–457. doi:10.1038/nsmb.2260. Trebec-Reynolds, D.P., et al., 2010. IL-1α and IL-1β have different effects on formation and activity of large osteoclasts. J. Cell. Biochem. 109, 975–982. doi:10. 1002/jcb.22476. Tu, S., et al., 2008. Overexpression of interleukin-1β induces gastric inflammation and cancer and mobilizes myeloid-derived suppressor cells in mice. Cancer Cell 14 (5), 408–419. doi:10.1016/j.ccr.2008.10.011. Varga, G., et al., 2016. Immune suppression via glucocorticoid-stimulated monocytes: a novel mechanism to cope with inflammation. J. Immunol. 193 (3), 1090–1099. doi:10.4049/jimmunol.1300891. Voronov, E., et al., 2014. The role IL-1 in tumor-mediated angiogenesis. Front. Physiol. 5 (14). doi:10.3389/fphys.2014.00114. Wang, D., et al., 2010. Structural insights into the assembly and activation of IL-1β with its receptors. Nat. Immunol. 11 (10), 905–911. doi:10.1038/ni.1925. Weber, A., et al., 2010. Interleukin-1 (IL-1) pathway. Sci. Signal. 3 (105), 1–6. doi:10. 1126/scisignal.3105cm1. West, N.R., et al., 2015. Emerging cytokine networks in colorectal cancer. Nat. Rev. Immunol. 15 (10), 615–629. doi:10.1038/nri3896. Zhuang, Z.-N., et al., 2016. IL-1 receptor antagonist inhibits pancreatic cancer growth by abrogating NF-kappaB activation. Clinic. Cancer Res. 22 (6), 1432–1444. doi:10.1158/1078- 0432.CCR- 14- 3382.