What determines the FDC polarization and GC size in affinity maturation

Germinal center (GC) is a particular biological structure produced for affinity maturation in the lymphoid follicle during the T-dependent immune response and is an important component of the humoral immune system. However, the impact of morphological features of the GC on antibody production is not clear. According to the latest biological experiments, we establish a spatiotemporal stochastic model to simulate the whole self-organization process of the GC including the appearance of two specific zones: the dark zone(DZ) and the light zone (LZ). We find that the development of light and dark zones in GC serves to maintain an effective competition among different cells and promote affinity maturation. On the other hand, by varying the GC size, a phase transition is discovered, which determines a critical GC volume for best performance in both the stochastic and the deterministic model. This critical volume is determined by the distance between the activated B Cell Receptor(BCR) and the target epitope of the antigen. The conclusion is confirmed in both the 2D and the 3D simulations and explains partly the variability in the GC size. Author summary Germinal center (GC) is an important component of the humoral immune system, which supports antibody affinity maturation and the generation of immunity memory. However, the impact of special morphological features of the GC on antibody production is not clear. According to the latest biological experiments, we establish a spatiotemporal stochastic model to simulate the whole self-organization process of the GC. We use the mixing index of different B cells to quantitatively describe the polarization in GC. With the increase of the mixing index, the affinity of plasma cells decreases gradually, even GC might collapse. Therefore, the development of light and dark zones in GC serves to maintain effective competition among different cells and promote affinity maturation. On the other hand, by varying the GC volume, a phase transition is discovered, which determines a critical GC volume for best performance in both the stochastic and the deterministic model. This critical volume is determined by the distance between the activated B Cell Receptor (BCR) and the target epitope of antigen. The conclusion is confirmed in both the 2D and the 3D simulations and explains partly the variability in the GC size.


Introduction
1 not only helps us understand the accelerated affinity maturation, but also suggests ways 46 in controlling the GC formation. 47 In this article, we first build a stochastic model based on which the impact of GC 48 morphology on affinity maturation is quantitatively studied. One interesting question is 49 why most groups of FDCs are polarized in GCs. From the simulation, we see that the 50 CXCL12-expressing reticular cells (CRCs) in DZ and FDCs in LZ together separate 51 somatic hypermutation(SHM) and mitosis from clonal selection, which serves to maitain 52 the competition of specific B cells and thus bestow efficient affinity maturation. Next, 53 we check the influence of the GC size on affinity maturation. Simulation results show 54 that the GC has an optimal volume for affinity maturation, determined by a phase 55 transition at a particular volume size beyond which the affinity of plasma cells is greatly 56 improved. In addition, this phase transition is first order, which is clearly seen in a 57 deterministic model derived from the stochastic one. More detailed study shows that 58 the phase transition point is determined by the distance in the gene space between the 59 activated B cell and the target epitope antigen in lymph nodes. The final basic question 60 is whether a 2D model can be used to represent the real GC reactions, so we extend the 61 model to 3D for comparison and find a positive answer.

63
Two models of the GC dynamics 64 In this paper, we establish a stochastic model as an extension of the one proposed in an 65 earlier article [16] being summarized in Appendix A. The simulation is carried out with 66 an accelerated version of the Gillespie algorithm [32]. More details are considered in the 67 current model together with some newly revealed regulatory biochemical reactions. For 68 example, we consider the entry and leaving of the Tfh cells in the GC and check the 69 dynamic equilibrium of the Tfh cells in the GC. Once fully activated B cells enter the 70 primary follicle, we deploy 100 activated T cells in paracortex, which ensures that the T 71 cells entering the GC accounts for 5% of the total number of cells in the GC. These 72 activated T cells actively respond to the CXCL13 secreted by FDCs [37], the 73 concentration of which is inversely proportional to the distance from the FDC. With 74 this attraction, these T cells gather around the primary follicle but are prevented from 75 entering by bystander B cells [29,30,38]. Thus, ICOSL expression by follicular 76 bystander B cells, engagement with ICOS signaling optimizes phosphoinositide-3 kinase 77 (PI3K)-dependent pseudopod dynamics, which promotes T cells persistent motility at 78 the T-B border and enables T cells to cross the border and enter the primary follicle. 79 According to a latest research from Professor Hai Qi, PD-1, another B7 family molecule, 80 though inhibitive, is highly expressed by T cells. They found that PD-1 inhibits the 81 entry of T cells into GC by suppressing PI3K , which was induced by the PD-L1 82 expressed on follicular bystander B cells [31]. Thus, in our model, we consider activities 83 of these co-inhibitory molecular pairs, leading to a dynamic equilibrium of the T cells. 84 The stochastic model 85 We use a rectangular grid to represent the entire GC, which has a spacing of 10µm 86 (Fig.1A), providing a platform for different cells to migrate randomly. We assume that 87 the same type of cells cannot occupy the same grid point, while different types of cells 88 can. In the model, CBs and CCs are regarded as the same type. The state of a cell is 89 described by a vector containing many elements, such as spatial location, protein 90 concentration, etc. During the GC evolution, each possible state change is regarded as a 91 reaction. When a reaction takes place, the corresponding vector element of the cell will 92 August 31, 2020 3/25 be updated. The signal integrator protein (SIP), a hypothetical protein that can 93 regulate the immune behavior of different B cells in the simulation (Appendix A). The 94 parameters and initial conditions of the model are shown in Table 2.

97
CCs could sense the concentration of CXCL13 secreted by FDCs. We assume that the 98 concentration of CXCL13 is inversely proportional to the distance from the FDCs, When T cells are activated in the paracortex, the generation rate of ICOS and PD-1 is assumed to be 0.02/h and 0.01/h. When T cells are in contact with bystander B cells, due to the activation of ICOS and PD-1, the speed of T cells changes according to: where n is the concentration of ICOS. m is the concentration of PD-1 in T cells. 105 P 1 = 4 × 10 −4 , A 1 = 8. Due to the physical blockage of bystander B cells, the rate at 106 which T cells enter a follicle is reduced. In our model, the speed of movement to the 107 follicle is reduced by 5µm/min. In the model, we use a four-dimensional lattice as the shape space to describe different 110 gene types of antibodies and the target epitope of the antigen ( Fig.2A). GC enriches the 111 BCR phenotype by rapid division and proliferation of B cells with high frequency 112 somatic mutations. In the stochastic model, the affinity of the antibody is assumed to 113 double or reduce by half with each mutation of B cells. The epitope of the antigen is 114 represented with a unique point -the origin φ * = 0 0 0 0 . A mutation is 115 represented by a jump to a nearest neighbor in the shape space. More explicitly, the 116 affinity of an antibody can be defined as: where φ is the coordinate of the BCR phenotype in the shape space, and the constant 118 k x = 0.0037, k y = 5.6, k z = 0.7. As can be seen from Fig.2B, the rate of taking up 119 antigens from FDCs is higher when the affinity is higher. With the highest affinity, the 120 CCs would take up one antigen in about every 4min.  The total number of activated Tfh cells 20 T The number of Tfh cells that are not contact with CCs - The number of CBs with affinity level i m The probability of no mutation 0.5 m 2 The probability of each mutation moves forward or move backwards 0.25 P i The total concentration of SIP with affinity level i - The number of CCs with affinity level i that have been rescued by T cells - The number of CCs with affinity level i that are binding with Tfh cells -R c The probability of CCs' recycling into the dark zone 0.8 The sum of the numbers of CCs and CBs with affinity level i - The SIP concentration of each B cell at level i - Based on the stochastic model described above, a series of equations can be written 123 down to describe the averaged dynamics in the GC. In the deterministic model, we set 124 the affinity of BCR to 10 levels (0 ∼ 9) from low to high, which determines the binding 125 rate of B cells to antigens. Correspondingly, the shape space is simplified a   138 where T 0 (= 20) is the total number of activated Tfh cells and T is the number of Tfh 139 cells that are not in contact with CCs. CC iT f h is the number of CCs with an affinity 140 level i and binding with Tfh cells.
which describe the binding of CCs with FDCs and the change of the antigen the degradation rate of ingested antigens. Ag i is the concentration of the ingested 153 antigens by CCs with affinity level i. R 6 is a rate defined in Appendix A6.
where F R i is the number of CCs with affinity level i that had contacted with FDC sites 155 and EN CT (=0.0001) is the possibility of CCs encountering with Tfh cells.

156
In the deterministic model, the positive feedback of ICOS-ICOSL and CD40L-CD40 157 and the impact of PD-L1-PD-1 molecular pairs are not considered. Although these 158 molecular pairs contribute to affinity maturation, the decisive factor is still the 159 concentrations of antigens absorbed by CCs. Therefore, this model is not identical with 160 the stochastic model and the impact of molecular pairs is not considered. Similar to what is in the stochastic model, the higher the antigen concentrations, the faster the 162 binding rate of CCs with Tfh cells. The reaction rates R 7 and R 8 : which describe the binding of CCs with Tfh cells and the help signals from Tfh cells. k T (= 0.05) is the dissociation rate of CCs and Tfh cells. F R iT is the number of CCs with affinity level i that are binding with Tfh cells. SP i is the concentration of help signals in CCs. k s (= 1) is the degradation rate of help signals in T R i .
which describe the process of CCs differentiating into CBs and into memory or plasma 164 cells. P B i represents output cells with different affinity levels, namely plasma cells and 165 memory cells.
which describes the concentration changes of SIP. k p is the degradation rate of total concentrations of SIP with affinity level i. In this model, the synthesis rate of SIP R 11 is affected by the concentration of help signals in CCs: where P is the concentration of help signals and Hill's coefficient m B = 4. The constant 167

169
We use the spatiotemporal stochastic model based on the Gillespie algorithm to 170 simulate the immune response in the GC within 21 days. When simulating in silico, we 171 can easily change the morphology of GC and analyze its impact on affinity 172 maturation [33,34]. We simulate the GC space using a 120*120 lattice with the grid 173 spacing being 10µm. The space was divided into two parts: the paracortex area and the 174 primary follicle area (Fig.1A).
At the beginning of the simulation, activated T cells and B cells moved to the T-B 176 boundary under the influence of chemokines. When these T cells meet an activated B 177 cell, this B cell will be activated fully and move to the primary follicle to generate GC. 178 At this time, FDCs will be polarized. Naïve B cells that exist in primary follicles make 179 up their boundaries as bystander B cells. With the proliferation and differentiation of B 180 cells (Fig.1B), LZ and DZ gradually appear in the GC, and the chemokines CXCL12 181 and CXCL13 play an important role in the whole process. The positioning of CBs in 182 the DZ depends on the expression of CXCL12 receptors, CXCR4 [37]. When CBs 183 differentiate into CCs, they express CXCR5, the receptor of CXCL13, which promotes 184 CCs motility towards the LZ. (Fig.1A, Fig.1B Six different radii, 50µm, 100µm, 150µm, 200µm, 300µm, 400µm, are used to 223 carry out the computation which produce different affinities for the output cells. From 224 Fig.4A, we can see from the simulation results that when the radius of GC is 200µm, it 225 is the most suitable for affinity maturation, and there is a phase transition in affinity 226 level when the radii is between 100µm and 150µm. When the phase transition occurs, 227 the standard deviation of the affinity reaches the maximum, which indicates that the 228 simulation results fluctuate greatly at this radius. With further increase of the radius, 229 the affinity increases, while the standard deviation decreases, and the simulation results 230 become more and more stable. When the radius exceeds the optimal radius, the affinity 231 starts to decrease gradually with the increase of the radius. The reason for this decrease 232 can be seen from Fig.4C. With the increase of the volume of GC, the accumulation of B 233 cell mutation becomes more and more, a large portion of the mutation is farther and 234 farther away from the target epitope of antigen, and the proportion of B cells with low 235 affinity increases, so it becomes more difficult to select CCs with high affinity. As shown in Fig.5, we simulate the development and affinity maturation of the GC 237 with the deterministic equation. On the 21st day, the average affinity of all output cells 238 was 0.37611. A similar phase transition can be seen in Fig.6C for the deterministic 239 model, where the volume of the GC is controled by the concentration of SIP in the 240 initial B cells. As a result, the initial SIP concentration could be used to represent the 241 GC volume. The higher the concentrations of SIP, the larger the volume of the GC 242 (Fig.6A). The jump of the affinity can also be observed when the concentration of SIP 243 in naive B cells is 251, which clearly indicates a phase transition (Fig.6C). When the 244 concentration of SIP is less than 251, the volume of GC will decay to zero with time, 245 while when the concentration of SIP is higher than 251, GC will reach a stable volume. 246 August 31, 2020 11/25 At 251, GC contained 110 B cells, which was the largest volume of GC (near the third 247 day). Similar results are obtained when the critical radius is 100µm in the stochastic 248 model (Fig.6D). The existence of the phase transition could also be explained 249 qualitatively in a straightforward manner. If the number of CBs generated by division 250 or proliferation is too small, there is not enough mutation accumulation to produce B 251 cells with high affinity. As a result, these B cells cannot get enough help from Tfh cells 252 and most of B cells will execute apoptosis, so the GC will eventually collapse. In order 253 to achieve affinity maturation, the GC needs a sufficient number of B cells to undergo 254 SHM. The phase transition point pins down the lower bound of this number in average. 255 Below, we will investigate what controls the critical number of B cells for the phase 256 transition. In the stochastic model, the initial B cells are assumed to have a genotype 257 1 1 1 1 , the distance between BCR phenotype and the target epitope of the 258 antigen is δ (φ, φ * ) = 4. This distance plays an important role in determining the 259 critical GC volume as will be seen next. If we set the initial BCR genotype to 260 1 1 0 0 and 2 2 2 2 respectively in the stochastic model, the distance to the 261 target epitope is δ (φ, φ * ) = 2 and δ (φ, φ * ) = 8. As can be seen from Fig.7A and B, the 262 affinities in these two cases do not reach the same level. When δ (φ, φ * ) = 2, even if the 263 volume of the GC is very small, the average affinity of output cells is greater than 0.5, 264 but when δ (φ, φ * ) = 8, the maximum affinity is 0.2, achieved at a radius as large as 265 400µm. Of course, by increasing the volume of the GC, it is possible to promote cell 266 affinity. With δ (φ, φ * ) = 2 increasing to δ (φ, φ * ) = 4, the GC radius with the highest 267 average affinity increases from 150µm to 200µm. Therefore, we may draw the 268 conclusion that the critical GC volume depends on the distance δ(φ, φ * ), which may be 269 checked experimentally.  When the volume of the GC is very small, it is important to find a mechanism to 271 ensure that plasma cells with high affinity could eventually be produced. In biological 272 experiments, it has been observed that GCs with different sizes grow in lymph nodes, 273 between which cells migrate now and then. Therefore, It is reasonable to speculate that 274 some B cells with high affinity may be able to move to smaller GCs to help them grow 275 and promote the production of high affinity plasma cells. To verify this hypothesis, in 276 both the deterministic and the stochastic model, we initially supply additional B cells  From Fig.8A, we can see that the affinity of the output plasma cells is 0.35, much 284 higher than 0.2, the affinity achieved without this additional input. However, if we 285 change the affinity of the ten B cells from 0.24 to 0.49 (Fig.8B), the affinity of plasma 286 cells is improved further to 0.5 and the standard deviation becomes much smaller. From 287 the results of Fig.6B and Fig.8, it is reasonable to say that B cells with higher affinity 288 that exit from GCs and enter other smaller GCs may help produce high affinity output 289 cells, and thus contribute to the collective maturation in the whole lymphoid. Simulation on the 3D lattice 291 All previous simulations are carried out on a two-dimensional lattice. However, realistic 292 biological experiments are performed and observed in three-dimensional space. 293 Therefore, whether the difference in dimensions would significantly change the results, 294 remains to be answered. First, we check the change of reaction rates going from 2D to 295 3D. In Meyer-Hermann's article [23], the values of the parameters in the 3D model are 296 directly computed from the 2D. The proliferation rate, differentiation rate, etc. must be 297 multiplied by 3/2 to ensure the results from both simulations are comparable, while the 298 concentrations of FDCs and T cells in GC should be kept constant. As a result, much 299 more FDCs and T cells (= 2700) are needed in three dimensions. However, the scaling 300 of the number of B cells is not so obvious and it is hard to make a direct comparison. So 301 a normalization of the B cell number is needed. As can be seen from Fig.9B, under the 302 same conditions, after the normalization that the number of B cells at each time should 303 be divided by the maximum number of B cells in the whole GC response, the B cells in 304 the 3D simulation seem to be greatly reduced in the later stage of the GC development. 305 It can be seen that the normalization of B cells in the GC in the 3D simulation is nearly 306 20% less than that in the 2D simulation under comparable conditions. The main reason 307 is that, in our simulation the number of B cells in the 3D is nearly 15 times of that in 308 2D. As a result even if the gene space and the mutation rate is the same, the mutation 309 range in the gene space of B cells expands. Many B cells stride randomly to a region 310 with very low affinity, so the affinity distribution of B cells changes. As shown in 311 Fig.9D, the proportion of B cells with affinity below 0.05 in 3D is significantly greater 312 than that in 2D. while 0.05 and 0.15, the portion in 2D is greater. As a result, many B 313 cells die quickly because their affinity is too low to get help from T cells. However, if we 314 continue to increase the number of T cells outside the GC from 2700 to 6000 in 3D, 315 arriving at a higher concentration than in 2D, B cells start to get more help from T 316 cells, resulting in a similar result to that of the 2D case. This also explains why large 317 volume of the GC in Fig.4A and Fig.4C leads to a decrease in the affinity.

318
In the 3D simulation, we also check how the radius of the GC alter the affinity 319 maturation. As can be seen from Fig.4B, a phase transition occurs at a radius of 50µm, 320 where the fluctuation rises significantly compared to the 2D case. With a radius of 321 35µm, the peak number of B cells in the GC is also about 110, being consistent with the 322 results from the 2D simulation (Fig.6D) and the deterministic model. Therefore, if the 323 number of the initial B cells is too small, the resulted few mutations are unable to 324 support a sufficient exploration of the gene space so that affinity maturation could not 325 be achieved before the decay of the GC. plays a very important role in this process. When FDCs are in an unpolarized state, the 332 GC will decay and disappear quickly, and the affinity level of the plasma and memory 333 cells will be very low. By slowly increasing the polarization of FDCs, the structure of 334 the light and the dark zone gradually appears, and the function of affinity maturation 335 gradually recovers, which does show the importance of polarization in the follicle. 336 However, in the process of simulation, because of cell crowdedness, we assume that cellular movement may sometimes be blocked. How much this simulation is consistent 338 with experimental observations remains to be explored. 339 Next, we explore how different volumes of the affect affinity maturation and whether 340 there is an optimal GC volume. By analyzing the average affinity of differentiated 341 plasma and memory cells in both stochastic and deterministic simulations, we find that 342 with the increase of volume, the average affinity of output cells undergoes a phase 343 transition beyond which a volume may exists that is most favorable for affinity 344 maturation. It turns out that the distance δ (φ, φ * ) in the gene space between the 345 activated B-cell and the target epitope of antigen determines the phase transition point 346 and thus the possible optimal volume of the GC. Conversely, different volumes observed 347 in the lymph node may be used as an indicator of this distance if the communication of 348 B cells between GCs could be excluded. Overall, if the number of cell divisions during 349 high frequency mutation is too small to produce enough strains with different affinities, 350 a devastating impact on the subsequent affinity maturation process will be resulted.

351
The affinity of the naive B cells is very low, being set at 0.06. By the end of the day, 352 the mean affinity of all generated plasma and memory cells rises to a value as high as  probability. In a two-dimensional plane, there are four directions to go. However, B cells cannot move to a lattice point occupied by another B cell. When FDCs are polarized in the follicle, the CCs will sense the concentration of CXCL13 secreted by FDCs. The CCs directed movement resulting from the chemoattractant allows the cells leaving the DZ before they die by apoptosis.
To favor the CCs with high affinity, we assume the following mechanism for competition. If one cell is in contact with the a FDC site, and its neighboring cell is not, these two B cells have chance to exchange their positions with probability proportional to the following quantity.
where ∆ A is the ratio of their BCR affinities. The constant c A = 2.0 gives the cell's basic competitiveness, and x 1 = 40.0, In this stochastic model, B cells have the opportunity to contact with T cells by presenting the pMHC molecule to the TCR. One B cell in contact with a Tfh cell may be replaced by a bystander B cell with more pMHC molecules on its surface. The exchange probability is : where, ∆ T is the concentration difference of pMHC and y 1 = 5.0, y 2 = 6.0, y 3 = 0.001.

A.2 Clonal proliferation of B cells
In the early stage of the immune response, B cells proliferate rapidly in the DZ, along with SHM. B cells with mutation have different BCR phenotype. The proliferation rate is approximately (9h) −1 and is controlled by SIP, expressed as: where n is concentration of SIP in the B cells, m 1 = 6, k 1 = 4.3519 × 10 −4 . At the initial stage, since the seeder B cells contain a high concentration of SIP (cl = 1000), they rapidly divide. SIP and other proteins are disseminated randomly to daughter cells after B cell proliferation. Over time, the cell division rate will decrease with the decrease of intracellular SIP concentration, so B cells need help signals from Tfh cells to synthesize more SIP.
In the simulation, the rate of B cell divisions depends not only on R 1 , but also on the number of vacancies in a neighborhood of the cell. Because each grid point can only be occupied by a cell of the same type, the number of B cells may seriously affect the rate of its division. For example, when there are eight B cells around a B cell, this B cell will not be able to divide if we assume that each one is fixed to the occupied lattice site. But this assumption is apparently not consistent with experimental observation. In fact, B cells can push other cells away for division. In order to cope with this problem, we assume that in a neighborhood with a radius of 100 µ m, if there are vacancies, the B cell may divide normally and the newly generated daughter-cell will be randomly put into to a vacancy in this neighborhood.

A.3 B cell differentiation
There is still no clear conclusion on how CBs differentiate into CCs. It has been stated based on experimental observation that the transition is independent of location. It August 31, 2020 21/25 might be triggered by LZ-derived signals such as the receipt of the T cell signal. In our model, the differentiation rate is controlled by the concentration of SIP. When the concentration of SIP in CBs was low, CBs began to differentiate into CCs. CCs, neither proliferating nor mutating express BCR on the surface of the cell membrane. They are in an activated state of apoptosis and have a life span of about 6 hours. In the model, the concentration of SIP in CBs determines the differentiation rate.
where n is the concentration of SIP in the B cells, m 2 = 5 is Hill's coefficient, k 2 = 0.0102. Proliferation of B cells leads to a decrease in the SIP concentration, which not only promotes the differentiate of CBs into CCs, but also increases the sensitivity of B cells to apoptosis. Therefore, the CCs with higher affinity are able to enter LZ to capture antigens and present them to T cells to reduce the rate of apoptosis. After the CCs receive the help signal from the Tfh cells, they start to synthesize and accumulate the SIP. These selected CCs may differentiate into CBs and return to the dark zone, where the mutation and selection is restarted. The differentiation rate R 3 of CCs to CBs is : where n is the concentration of the SIP, m 3 = 4 is the Hill's coefficient, k 3 = 0.0011.

A.4 Apoptosis of cells
According to experiments in vivo, CCs can survive in the LZ for about 6 ∼ 16 hours. Because SHM makes BCR affinity differ from one another, BCR goes through the selection of antigens presenting on FDCs and only those B cells expressing high affinity antibodies can effectively seize antigens from FDCs. Those with low affinity will be eliminated due to their high sensitivity to apoptosis. It is worth noting that the B cell apoptotic rate is essential in the regulation of GC during the immune response. If the rate of apoptosis of B cells is too high, B cells will all die in a short time. A high immunity level is not achievable. Conversely, if the rate of apoptosis of B cells is too slow, B cells will occupy most lattice points in the GC. As a result, the SHM cannot proceed smoothly, and finally slow down of BCR affinity maturation. During the simulation, the state vector of apoptotic B cells will be directly deleted and removed from the computer memory. In the model, we assume that the apoptotic rate increases with the decrease of SIP, so the rate R 4 is : where n is the concentration of the SIP, m 4 = 5 is the Hill's coefficient, k 4 = 0.05.

A.5 Interaction between CCs and FDCs
Follicular dendritic cells (FDCs), the most powerful antigen-presenting cells (APCs) in GC, efficiently absorb, process, and present antigens. The presented antigen binds to BCR to form an immunological synapse and plays an important role in activating B cells. B cells with high affinity will endocytosis antigens, and then within the B cell, antigens will be decomposed into polypeptide fragments, part of which then bind to the MHC class II molecule to form a complex and are transferred to the B cell surface. These complexes can be recognized by Tfh cells. In the model of this paper, CCs enter the lZ of GC, and may undergo three states: "not selected", "contact with antigen" and "antigen selection":

A.5.1 Unselected
CCs that are not in contact with antigens on FDCs will not be recognized by Tfh cells. Therefore, the concentration of SIP in B cells do not increase. The apoptosis was highly sensitive to lower concentration of SIP in B cells, so CCs will be eliminated by apoptosis at a rate of R 4 .

A.5.2 Contact with antigens
The "contact with antigens" CCs indicate that BCR has successfully contacted with antigens in FDCs, but there is no endocytosis antibody -antigen immune complex. We assume that the binding rate is related to the affinity of the BCR, and the CCs remain bound to the antigens for about 5min. The antigen binding rate R 5 is recorded as: where k a = 2.696, k b = 0.8, k c = 0.25. δ (φ, φ * ) represent the distance between the BCR and the antigen in the gene space δ (φ, φ * ) = 4 i=1 |φ i − φ * i |.

A.5.3 Selected by antigens
CCs selected by antigens indicate that BCR successfully ingested antigens from FDCs. In this model, the reaction rate of B cells to take up antigens is related to the affinity of BCR. It is given as: where k x = 0.0056, k y = 5.6, k z = 0.7. With the highest affinity, a CC will take up one Ag about every 4min. As you can see from Fig.2A, the rate of uptaking Ag from FDCs is higher when the affinity is higher.

A.6 Interaction between CCs and Tfh cells
Professor Hai Qi used TIRF microscopy to observe Tfh cells entangled with B cells.
During T-B entanglement, two co-stimulatory molecules of ICOS-ICOSL and CD40L-CD40 constitute a positive feedback network for improving affinity maturation of BCR. Recently, they found the role of another pair PD-1-PD-L1, the programmed cell death-1(PD-1) expressed by Tfh cells that controls the location and function of Tfh cells. It may dampen TCR signaling and thereby reduce the ligand sensitivity of Tfh cells, which should increase the overall stringency of selection in GCs. In real biological experiments, PD-L1 deficient B cells can get more help from Tfh cells. Based on these data, it was shown that the interaction of PD-1 and PD-L1 could also increase the affinity of BCR. These findings should be reflected in the stochastic model. Similarly, during the simulation, CCs after uptaking antigens could dwell in three different states: "not selected by T cells", "contacted by T cells" and "selected by T cells": A.6.1 Not selected by T cells CCs that are "not selected by T cells" indicate those that have not obtained help signal from T cells since entering the LZ. Therefore, these CCs will be difficult to survive in the LZ unless they contact with a T cell during the wandering in the LZ.