|Home | About | Journals | Submit | Contact Us | Français|
The formation of the protein segregation structure known as the “immunological synapse” in the contact region between B cells and antigen presenting cells appears to precede antigen (Ag) uptake by B cells. The mature B cell synapse consists of a central cluster of B cell receptor/Antigen (BCR/Ag) complexes surrounded by a ring of LFA-1/ICAM-1 complexes. In this study, we used an in silico model to investigate whether cytoskeletally driven transport of molecules toward the center of the contact zone is a potential mechanism of immunological synapse formation in B cells. We modeled directed transport by the cytoskeleton in an effective manner, by biasing the diffusion of molecules toward the center of the contact zone. Our results clearly show that biased diffusion of BCR/Ag complexes on the B cell surface is sufficient to produce patterns similar to experimentally observed immunological synapses. This is true even in the presence of significant membrane deformation as a result of receptor–ligand binding, which in previous work we showed had a detrimental effect on synapse formation at high antigen affinity values. Comparison of our model’s results to those of experiments shows that our model produces synapses for realistic length, time, and affinity scales. Our results also show that strong biased diffusion of free molecules has a negative effect on synapse formation by excluding BCR/Ag complexes from the center of the contact zone. However, synapses may still form provided the bias in diffusion of free molecules is an order-of-magnitude weaker than that of BCR/Ag complexes. We also show how diffusion trajectories obtained from single-molecule tracking experiments can generate insight into the mechanism of synapse formation.
The “immunological synapse” is a membrane–protein segregation structure that forms during contact between a lymphocyte and an antigen presenting cell (APC) during recognition of antigen by the lymphocyte. The function and formation mechanisms of the immunological synapse are among the least understood aspects of lymphocyte activation. Originally observed to form during antigen recognition by T cells,17,20,24,35 immunological synapses have also been observed during contact between B cells and APCs.2 During synapse formation in B cells, the B cell receptors (BCR) bind antigen (Ag) on the APC surface and the resultant BCR/Ag complexes cluster at the center of the contact zone, while the integrin Lymphocyte function-associated antigen-1 (LFA-1) on the B cell surface binds its ligand Intercellular adhesion molecule-1 (ICAM-1). The LFA-1/ICAM-1 complexes surround the central cluster of BCR/Ag complexes, resulting in a concentric pattern.2,6,8 It is widely believed that the immunological synapse modulates intracellular signaling and immune cell response.7,11,15,22
A considerable modeling effort has been expended to understand the formation mechanisms of the immunological synapse, particularly for the T cell synapse,4,9,10,16,22,23,25,26,33 and to a lesser extent the B cell synapse,15,18,31,32 although many aspects of synapse formation still remain unresolved. In T cells, differences in equilibrium bond length between TCR/MHCp and LFA-1/ICAM-1 complexes are thought to be sufficient to generate immunological synapses (the so-called “topographic model”).25,26 However, experimental studies involving T cells show that receptor accumulation at the T cell-APC interface is also driven by the cell’s cytoskeleton.5,13,19,34
Even though immunological synapse patterns observed in B cells resemble the canonical T cell synapse pattern, the significantly larger size of the BCR (~25 nm1) compared to the TCR/MHCp bond (~15 nm4) means that differences in equilibrium bond length between BCR/Antigen and LFA-1/ICAM-1 complexes (equilibrium bond length ~40 nm4) are less likely to account for synapse formation in B cells. In addition, B cells recognize antigen over a considerably wider range of affinity values (KA = 106–1011 M−1) than do T cells (KA = 106–108 M−1). Other significant differences between B cells and T cells include receptor valency (the BCR is bivalent, as compared to the TCR, which is monovalent), and that the density of B cell receptors on the B cell surface is at least an order of magnitude greater than that of T cell receptors on the T cell surface. Lastly, the number of experimental studies that focuses on exploring the molecular basis of B cell synapse formation is small compared to such studies for T cells, making it even more difficult to reach a consensus on the formation mechanism of immunological synapses in B cells.
We have previously studied B cell synapse formation using Monte Carlo simulations.31,32 Our modeling work indicated that it is difficult for high affinity antigens (KA ≥ 108 M−1) to cluster into a synapse by purely passive mechanisms, such as a difference in bond length between the BCR/Ag and LFA-1/ICAM-1 complexes, in the occurrence of significant membrane deformation.31 We also showed that BCR molecules undergo sub-diffusive motion upon binding antigen, making it difficult for them to cluster in a synapse pattern within biologically realistic time-scales.32 It thus seems likely that synapse formation in B cells is driven by an active, signaling-driven mechanism that involves the actin cytoskeleton.28,29 A recent study of B cell activation also showed that BCR accumulation at the synapse is compromised when the cytoplasmic signaling domain of BCR molecules has been truncated so as to inhibit their ability to signal and attach to the actin cytoskeleton.15
In this study we use a Monte Carlo simulation procedure to investigate whether cytoskeletally driven transport of receptor molecules is a potential mechanism of synapse formation in B cells. Cytoskeletally driven transport toward the center of the contact zone is simulated in a computationally efficient manner, by biasing the random diffusion of membrane bound molecules toward the center of the contact zone. In a series of in silico parametric experiments, we vary the strength of the bias in diffusion of the various membrane-bound species so as to generate both qualitative and quantitative insight into the nature of the synapse formation mechanism.
Our results indicate that a bias in the diffusion of BCR/Ag complexes toward the center of the cell–cell contact zone is a sufficient mechanism of synapse formation in situations where mechanisms that depend on differences in bond properties between BCR/Ag and LFA-1/ICAM-1 fail to generate synapses. However, our results also show that a strong bias in diffusion of LFA-1/ICAM-1 complexes and free molecules has a detrimental effect on synapse formation, although synapse formation can still occur if the bias in diffusion of these species is weak compared to that of BCR/Ag complexes. We also show how molecular diffusion trajectories obtained from single-molecule tracking experiments can be used to test for the presence of directed transport. If such transport is present, we show that the mean distance of receptor–ligand complexes from the center of the contact zone can be used to determine precisely which species are affected by the transport mechanism, as well as the relative strength of directed transport of the various species. Specifically, we show that if cytoskeletally driven transport affects BCR/Ag complexes more strongly than LFA-1/ICAM-1 complexes, or free receptor molecules, the respective mean distance from the center differs by an order of magnitude. Recently, Tolar et al.28,29 have studied receptor diffusion using single-molecule tracking during the course of B cell synapse formation. Our results can thus be used synergistically with parallel single molecule tracking experiments such as those in Tolar et al.28,29 to elucidate the mechanism of synapse formation in B cells.
We use a Monte Carlo procedure similar to our previous work.31,32 Receptor and ligand molecules are randomly sampled to undergo reaction or diffusion according to specific probabilities. A distinguishing feature of our method is a mapping between the probabilistic parameters of the Monte Carlo simulation and their physical counterparts. This makes it possible to compare our model’s results to those of physical experiments to within an order of magnitude.
We model a B cell–lipid bilayer system similar to that used in B cell synapse formation experiments,8,15 where the lipid bilayer mimics an antigen presenting cell surface. The simulation domain consists of a 3 × 3 μm2 square area on the bilayer surface and its vertical projection on the B cell surface, as shown in Fig. 1. This area is chosen so that the vertical separation distance at the edges of the domain (~400 nm) is an order of magnitude greater than the typical receptor–ligand bond (~40 nm), ensuring that the simulation domain encompasses the entire area where receptor–ligand binding is possible. This also allows us to assume a zero net flux boundary condition, which in our model is simulated by reflective boundaries. The lipid bilayer and B cell surface are modeled as discrete, N × N Cartesian lattices. We assume the B cell is initially spherical in shape (so as to minimize free energy), and the vertical separation distance z between the two surfaces at any given point (x, y) is thus given by
Only one molecule can occupy a node, so we choose a nodal spacing of 10 nm, which is approximately equal to a membrane protein’s exclusion radius (resulting in N = 300 nodes). The exception are BCR molecules, which being bivalent, have a width of ~25 nm1 and thus occupy three nodes, with either a horizontal or vertical orientation on the lattice. For the radius of B lymphocytes we use RB = 6 μm and z0 = 40 nm.
At the start of a simulation run, molecules are uniformly distributed over the two surfaces at random. The molecular species simulated are BCR and LFA-1 on the B cell surface, and their ligands, antigen and ICAM-1, on the bilayer surface. At each time step in the simulation, the molecules are individually sampled at random to undergo either diffusion or reaction events, determined by means of a coin toss with probability 0.5.
If a molecule has been selected to undergo reaction, we first check the facing node on the opposite surface for a binding partner. If that is the case, a random number trial with probability pon(i) is performed to determine if the two molecules will form a receptor–ligand complex. BCR molecules are able to bind two antigen molecules, one on each end node (but not the middle node). Thus, if a BCR molecule is selected for a reaction, an additional coin toss is performed to pick one of the end nodes, and the bilayer surface opposite the chosen node is checked for a free antigen molecule. Sometimes a BCR molecule may have bound an antigen molecule on one Fab domain and have the other Fab domain free, forming a BCR/Ag complex. If the free Fab domain is selected, the reaction proceeds as described above, which may result in a second antigen molecule binding to the BCR/Ag complex (forming a BCR/Ag2 complex). If the Fab domain with the bound antigen is selected, the BCR/Ag complex may dissociate into its component molecules with probability poff(i). Three reversible reactions are thus possible: LFA-1 + ICAM-1 ↔ LFA-1/ICAM-1, BCR + Ag ↔ BCR/Ag, and BCR/Ag + Ag ↔ BCR/Ag2. The binding and dissociation probabilities for the two reactions involving antigen are assumed to be the same and thus the subscript i refers to the BCR/Ag reactions when i = BA and the LFA-1/ICAM-1 reaction when i = LI. The overall sampling rate for reaction or dissociation is thus 0.5 × pon(i) or 0.5 × poff(i).
We assume the probability of bond formation depends on the vertical separation distance z in accordance with the linear spring model.3,12 Replacing the rate constant kon with the probability pon, we obtain the following probability density function
The bond is modeled as a mechanical spring with stiffness κ′ and equilibrium length zeq, while kB denotes the Boltzmann constant (1.34 × 10−23 J/K) and T the temperature (~300 K). Similarly, the dissociation probability of a receptor–ligand complex, poff(i), is given by
Without loss of generality, we can set κ(i) = 2κ′(i) so that the exponential in Eq. (3) is the same as that in Eq. (2) but with a positive sign in front. In contrast to pon(i), poff(i) is a minimum at z = zeq.
The quantity PA(i)(z) defined in Eq. (4) is analogous to the overall receptor–ligand affinity, and consists of both the intrinsic affinity PA(i)max and the bond stiffness κi. Varying ponmax and poffmin while keeping the ratio PA(i)max constant changes the time scale of the simulation, but not the equilibrium behavior.
If a molecule has been selected to undergo diffusion, a random number trial with probability pdiff(i) is used to determine whether the diffusion move will occur successfully. The overall sampling rate for diffusion is thus 0.5 × pdiff(i). Although the probability of diffusion can be different for all seven species present in the simulation (free BCR, free antigen, free LFA-1, free ICAM-1, BCR/Ag complexes, BCR/Ag2 complexes, LFA-1/ICAM-1 complexes), we assume for simplicity that free molecules in the bilayer (antigen, ICAM-1) diffuse with the same probability (denoted as pdiff(FB)), molecules on the cell surface diffuse with probability pdiff(FC), and receptor–ligand complexes with probability pdiff(C).
If the trial with probability pdiff(i) is successful, one of the four neighboring nodes is selected at random for the molecule to diffuse to. Because two molecules are not allowed to occupy the same node, the molecule will only move if the target node is unoccupied. For the case of BCR molecules, three nodes need to be free for the molecule to diffuse in the direction transverse to its length, while only one free node is needed in order for it to diffuse along its length. In the case of complexes, the target nodes on both surfaces need to be free (two nodes for monomeric LFA-1/ICAM-1 complexes, two or four for BCR/Ag complexes, and three or five for BCR/Ag2 complexes).
Because of the intricacies and computational cost involved in explicitly modeling the cytoskeleton, cytoskeletally driven motion toward the center of the contact zone is simulated in a computationally efficient, indirect manner by biasing the diffusion of molecules toward the center (see Appendix for details). We define a biasing factor η that we multiply pdiff(i) by if the target node is closer to the center of the simulation domain than the molecule’s current location. The case η = 1 corresponds to purely random motion, with biased motion toward the center increasing with larger η.
Our model includes changes in the originally spherical cell shape as a result of receptor–ligand binding. This is modeled by changes in the vertical separation distance z according to the approach used in Qi et al.,25 Raychaudhuri et al.,26 Chakraborty et al.,9 and Tsourkas et al.,31 which is given by the following equation
The quantities CBA and CLI refer to the local concentrations of BCR/Ag and LFA-1/ICAM-1, respectively, while γ and β represent membrane tension and bending rigidity, respectively. The constant M relates the time scale of membrane deformation to that of receptor–ligand binding, such that for small M, the membrane will essentially retain its shape for the duration of the simulation. Because the length scale of membrane deformation (set by (β/γ)1/2) is considerably larger than that of a protein’s exclusion radius (~100 nm instead of ~10 nm), for the purpose of calculating z we coarse-grain the N × N membrane surface lattice into 10 node × 10 node subdomains over which z is constant. At the end of each time step, the values of CBA and CLI in each subdomain are calculated and entered into a discrete version of Eq. (5), which is then solved numerically.
The molecular population is sampled S times for diffusion or reaction during every time step. The number of trials S is set equal to the total number of molecules (free and complex) present in the system at the beginning of each time step, and the simulation is run for a number of time steps T. The vertical separation z between the surfaces is updated at the end of each time step according to Eq. (5). A schematic of our Monte Carlo algorithm is shown in Fig. 2.
The parameters used in our model are listed in Table 1. Parameter values found in the literature are given on the left side of Table 1, while the appropriately mapped forms used in our model are listed on the right side of Table 1.31 Parameters whose values vary during experiments (such as BCR/Ag affinity and antigen concentration) are also varied in our simulations. We were not able to find measured values for the stiffness and extracellular domain length of the BCR/Ag bond, κBA and zeq(BA), and also for the membrane deformation time scale parameter M. We therefore assumed the stiffness and extracellular domain length of the BCR/Ag bond to be roughly similar to the corresponding values for the LFA-1/ICAM-1 pair, while we set M to a value of 10−12 m4/J s, which is the minimum value for which significant membrane deformation will occur within the time scale of synapse formation.
The diffusion coefficient of free receptor molecules in a cell membrane is in the range of ~0.01–0.1 μm2/s,14 while we estimate that it is an order of magnitude greater for molecules in a lipid bilayer. Since free antigen and ICAM-1 on the bilayer are the fastest diffusing species, we set pdiff(FB) = 1 and accordingly pdiff(FC) = 0.1. Due to the difficulty of experimentally measuring the diffusion coefficient of receptor–ligand complexes, we have not been able to find measured values for this quantity. B cell synapse experiments, however, show that a significant loss in mobility occurs upon antigen binding.28,29 Consequently, we use a value of pdiff(C) = 0.01 in our simulations. In this case, the bias in diffusion will affect free molecules and receptor–ligand complexes differently, with slow diffusing complexes being proportionately more affected (see Appendix). We carry out in silico experiments in which we vary the ηi for free BCR, BCR/Ag complexes, free LFA-1, and LFA-1/ICAM-1 complexes. Free antigen and ICAM-1 molecules, here modeled as being on artificial lipid bilayer, are always assumed to have ηi = 1.
Because some of the parameters of our model are probabilistic in nature and therefore dimensionless, it is necessary to map them onto physical quantities in order to physically interpret the results. Two such mappings are necessary: One which maps the probabilistic affinity PAmax to the association constant KA and one which maps the size of our model’s time step to physical time by relating pdiff to the physical diffusion coefficient D.
To map values of PAmax onto corresponding values of KA, we make use of the fact that at kinetic equilibrium, the two-dimensional association constant, KA(2D), can be obtained from the following relation3,4,21
Here C refers to the concentration (molecules/area), Ncomplex is the number of complexes formed at equilibrium, while Nfree(1) and Nfree(2) refer to the number of free molecules present at equilibrium. To map PAmax to KA(2D), we run our simulation for a particular value of PAmax to obtain Ncomplex, Nfree(1), and Nfree(2), and calculate KA(2D) from Eq. (6). Because KA is usually given in units of M−1, we multiply KA(2D) by the effective confinement length in the manner of Bell,3 for which we use the thickness of a cell membrane (~10 nm). When this is done, we obtain the following linear relationship between PAmax and KA(3D)31
Next, we establish the mapping of our model’s time scale to physical time by mapping the probability of diffusion pdiff to the diffusion coefficient D by means of direct simulation.31 From these simulations, pdiff = 1 maps to Dsim = 1 (nodal spacing)2/time step, or 10−4 μm2/time step for a 10 nm nodal spacing. By equating Dsim = 10−4 μm2/time step to the diffusion coefficient of the fastest diffusing species (since pdiff cannot exceed 1), in this case free antigen and ICAM-1 on the bilayer (~1 μm2/s), we calculate that the size of our simulation’s time step is 10−4 s, from the following relation
In Fig. 3, we plot the distribution of bound antigen (green) and bound ICAM-1 molecules (red) after T = 106 time steps (100 s) for BCR/Ag affinity KA(BCR/Ag) = 1010 M−1 (highest physiological value) and A0 = 100 antigen molecules. BCR/Ag complexes are scattered uniformly throughout the zone of contact when there is no bias in diffusion toward the center (ηBA = 1), with the result being no synapse formation in this case (Figs. 3a and and3e).3e). This is to be expected in the absence of driving factors for BCR/Ag complex clustering, such as a difference in bond length or bond stiffness between BCR/Ag complexes and LFA-1/ICAM-1 complexes.4,18,25,31 A significant increase in BCR/Ag complex clustering is observed when their biasing factor ηBA is set to 1.1 (Figs. 3b and and3f),3f), which increases further as the biasing factor ηBA is increased to 1.2 and 1.5 (Figs. 3c, c,3g3g and and3d,3d, d,3h,3h, respectively). The patterns formed with ηBA = 1.2 and 1.5, in particular, are highly reminiscent of experimental B cell synapses observed at similar values of BCR/Ag affinity and antigen concentration.2,6,8,15 These results therefore indicate that biased diffusion of BCR/Ag complexes is sufficient for the generation of B cell synapses, even when significant membrane deformation occurs. It is interesting to note that even a modest bias in diffusion of BCR/Ag complexes toward the center is sufficient to induce substantial clustering and generate good quality synapses.
Biased diffusion of BCR/Ag complexes also results in synapse formation for the entire physiological range of BCR/Ag affinity (KA = 106–1010 M−1) as well as antigen concentration (A0 = 10–100 molecules/μm2). However, our simulation failed to generate synapses below a BCR/Ag affinity value of KA = 106 M−1. In Figs. 4a and and4c,4c, a synapse forms when BCR/Ag affinity is set to KA = 106 M−1 with A0 = 1000 antigen molecules (~100 molecules/μm2), while in Figs. 4b and and4d,4d, where BCR/Ag affinity is set to KA = 10−5 M−1, we see that no synapse forms, as the number of BCR/Ag complexes formed is too low. This result persists even if we increase the concentration of antigen molecules and the strength of the bias in diffusion of BCR/Ag complexes. These results are particularly interesting in light of the fact that experimental investigations of synapse formation in B cells show that synapses do not form below a threshold BCR/Ag affinity of KA = 106 M−1, regardless of antigen concentration.2, 6, 8, 15
We also investigate the possibility that other membrane-bound species, such as free BCR molecules, free LFA-1 molecules, and LFA-1/ICAM-1 complexes, are transported toward the center of the contact zone by the cytoskeleton. In Fig. 5, the diffusion bias of BCR/Ag complexes is fixed at ηBA = 1.1, while the diffusion bias of free BCR, free LFA-1, and LFA-1/ICAM-1 complexes (ηB, ηL, ηLI, respectively), is varied as shown in the figure. When the diffusion of free BCR, free LFA-1 and LFA-1/ICAM-1 complexes is weakly biased (ηB = ηL = ηLI = 1.01, 1.02, Figs. 5a and and5b,5b, respectively), synapse formation is observed. However, if the bias in diffusion of these three species is increased further, BCR/Ag molecules fail to cluster and no synapse is observed (Figs. 5c and and5d).5d). This is due to crowding at the center of the contact zone by these species, in particular the free molecules. A strong bias in the diffusion of free molecules will lead to rapid crowding of the central zone by those molecules before significant BCR/Ag binding has had time to occur, and thus exclusion of BCR/Ag complexes from the center. It should be noted that a weaker bias in the diffusion of free BCR molecules compared to BCR/Ag complexes is probable, especially if there is an order-of-magnitude difference between the diffusion coefficient of free BCR molecules and BCR/Ag complexes, which is the case here.31,32 LFA-1/ICAM-1 complexes and free LFA-1 molecules may also have a weaker bias in diffusion compared to BCR/Ag complexes, possibly due to an intrinsically lower probability of attachment to the cytoskeleton.
We believe that computational modeling of B cell synapse formation can aid in devising experimental methods to test for the presence of cytoskeletally driven transport of membrane-bound molecules toward the center of the cell–cell contact region. In Fig. 6, we plot trajectories of individual antigen molecules during two sample runs of our simulation. Figure 6a corresponds to the case of Fig. 3a, where the diffusion of all membrane-bound species is purely random and no bias in diffusion is present (η = 1 for all species), while Fig. 6b corresponds to the case of Fig. 5a, where the diffusion bias of BCR/Ag complexes is an order of magnitude greater than that of the other three membrane bound species (ηBA = 1.1, ηB = ηL = ηLI = 1.01). Black stars in the figure indicate the location of the tracked molecule at the beginning of the simulation, while blue stars indicate the position at the end of the simulation. The black circles in the figure indicates the area where the vertical separation distance z is small enough for receptor-ligand binding to occur. Antigen molecules that are outside of the circles at the end of the simulation can be assumed to be free. As can be seen in the figure, there is a clear qualitative difference between the trajectories in Figs. 6a and and6b,6b, with the antigen molecules in Fig. 6b distributed considerably closer to the center than those in Fig. 6a. By performing a controlled experiment using wild-type and signaling-defective BCR, as in Fleire et al.15 we believe it would be possible to test for cytoskeletally driven transport of membrane-bound molecules by tracking individual antigen molecules and observing their trajectories.
In Fig. 7 we plot the mean distance from the center of BCR/Ag complexes, free BCR molecules, free LFA-1 molecules, and LFA-1/ICAM-1 complexes at the end of simulation runs corresponding to Figs. 3 and and55 (Figs. 7a and and7b7b respectively). In Fig. 7a, where only the diffusion of BCR/Ag complexes is biased, the mean distance from the center of BCR/Ag complexes drops by a factor of 2 when the biasing factor ηBA is increased from a value of 1.0 to a value of 1.1. The mean distance from the center of LFA-1/ICAM-1 complexes, free LFA-1, and free BCR remains constant as the diffusion of these species is not biased in this case. The mean distance from the center of free BCR and LFA-1 of the order of ~1.1–1.2 μm corresponds to the average distance from the center for uniformly distributed molecules (a circle of this radius covers half of the 3 × 3 μm2 contact area). From Figs. 7a and and3,3, we can state that a synapse is observed when the mean distance from the center of BCR/Ag complexes is about half that of LFA-1/ICAM-1 complexes (0.2 μm compared to 0.4 μm). We also note that the mean distance from the center of free BCR molecules is approximately an order of magnitude greater than that of antigen molecules (for the parameter values of Fig. 3, all antigen molecules are bound to BCR at the end of the simulation, thus the mean distance from the center of the antigen molecules is equivalent to that of the BCR/Ag complexes).
In Fig. 7b, the mean distance from the center of BCR/Ag complexes increases as the diffusion bias factor of free BCR, free LFA-1, and LFA-1/ICAM-1 complexes increases, and the mean distance from the center of free BCR and free LFA-1 molecules decreases. This is in spite of the fact that the diffusion bias factor of BCR/Ag complexes is fixed at ηBA = 1.1. The mean distance from the center of LFA-1/ICAM-1 complexes decreases, albeit weakly, as the effect of increasing diffusion bias (ηLI) is offset by increasing crowding at the center of the contact zone by free BCR and LFA-1 molecules. We also note that when the diffusion bias of free BCR and LFA-1 is strong, the mean distance from the center of these molecules is comparable to that of BCR/Ag complexes, while when the bias of free BCR and LFA-1 molecules is weak, the mean distance from the center is approximately an order of magnitude greater than that of BCR/Ag complexes. These data suggest that it is possible to obtain insight into the nature of the cytoskeletal transport mechanism by calculating the mean distance from the center of the contact region of the various species involved in B cell synapse formation.
In this study, we used a kinetic Monte Carlo simulation method to investigate whether directed transport of molecules toward the center of the B cell/bilayer interface is a potential mechanism of immunological synapse formation in B cells. A distinguishing feature of our method is the development of a mapping between probabilistic parameters of the Monte Carlo simulation and their physical counterparts, thereby allowing quantitative comparison of our model’s results to those of biological experiments.27,31,32 Significantly, the length and time scales of synapse formation in our model matched experimentally observed length and time scales of B cell synapse formation. The synapses formed in our model are of the order of ~1–1.5 μm in diameter, which is comparable to the diameter of physiological B cell synapses, while the formation time of ~100 s predicted by our model matches relatively well with the experimentally observed synapse formation time of 1–2 min. In addition, our model also reproduced the existence of a threshold value of BCR/Ag affinity for synapse formation of KA = 106 M−1. As in synapse formation experiments,6,8,15 our model did not produce synapses below this critical value of BCR/Ag affinity, thereby matching the experimentally observed affinity range of synapse formation.
We modeled directed transport of receptors in an implicit manner, by biasing the diffusion of molecules toward the center of the B cell/bilayer interface. This approach, rather than an explicit simulation of molecular attachment to the cytoskeleton, was chosen based on computational efficiency considerations. In the Appendix, we show that biased diffusion is mathematically equivalent to explicit simulation of attachment and detachment to the cytoskeleton. Our results reveal that biased diffusion of BCR/Ag complexes is sufficient to produce synapse patterns similar to those observed in biological experiments, even for high affinity antigens. In the absence of such a transport mechanism, BCR molecules show sub-diffusive behavior in the early phase of synapse formation, which inhibits rapid clustering of BCR/Ag complexes, especially for high affinity antigens.32 In the late phase of synapse formation, significant membrane deformation makes it difficult for synapses to form for high affinity antigens in the absence of a transport mechanism.31 These results are thus all the more important in light of the fact that BCR can bind antigens with high affinity (up to 1010 M−1), and that there is experimental evidence to suggest that significant membrane deformation does indeed occur during B cell synapse formation.2,6,8,15
Our results lie in contrast to theoretical studies of T cell synapse formation where it was shown that a difference in equilibrium bond length between TCR/MHCp and LFA1/ICAM1 complexes was enough to segregate the two receptor–ligand pairs into an immunological synapse pattern.25,26 However, it has been shown during T-cell synapse formation experiments that the time needed for an antibody-coated bead to traverse half a cell circumference is approximately 6 min (for a ~10 μm diameter T cell).34 Assuming a diffusion coefficient of the order of ~0.1 μm2/s for receptors on a cell surface,14 a simple calculation (time~distance2/diffusion coefficient) shows that the time needed for receptor to traverse half a cell circumference by pure diffusion is much larger than the time indicated in the work of Wülfing and Davis.34 Thus, it seems likely that active transport of receptors assists the bond length difference-mediated synapse formation mechanism in T cells. In B cells, however, the absence of a difference in equilibrium bond length between BCR/Ag and LFA-1/ICAM1 makes these species’ spontaneous segregation into a synapse pattern difficult. Moreover, in previous work, we also showed that canonical synapse formation by purely passive mechanisms based on differences in the properties of the BCR/Ag and LFA-1/ICAM-1 bonds (e.g., equilibrium bond length, bond stiffness, affinity) is not possible, especially for high affinity antigens on a cell membrane that undergoes deformation.31 However, the results of this study show that directed transport of BCR/Ag complexes to the center of the cell–cell interface is capable of generating canonical synapse patterns even in situations where mechanisms based on differences in bond properties between BCR/Ag and LFA-1/ICAM-1 fail to do so. Importantly, recent B cell synapse formation experiments show impaired synapse formation in signaling-deficient B cells.15 The likely explanation for this is that signaling-deficient B cell receptors cannot generate the signal needed for cytoskeleton-mediated transport of B cell receptors to the synapse.
Our results also indicate that biased diffusion of the other membrane-bound species involved in B cell synapse formation (free BCR, free LFA-1, LFA-1/ICAM-1 complexes) hinders synapse formation unless the bias in the diffusion of BCR/Antigen complexes is significantly stronger than that of the other species. A similar, differential transport mechanism has been shown to be effective in experimental studies of T cell synapse formation.19
We show that it is possible to gain insight into the mechanism of synapse formation by obtaining the trajectories of individual antigen molecules over the course of synapse formation, as well as by calculating the mean distance from the center of the various species. The trajectories can potentially reveal the presence or absence of cytoskeletal transport, while the mean distance from the center can reveal whether cytoskeletal transport affects the various species differently. We believe that an iterative process of computational modeling and physical experimentation can lead to a full understanding of the mechanism of immunological synapse formation in B cells.
PT and SR are supported from NIH Grant AI074022.
Open Access This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
Due to the computational cost of explicitly modeling the B-cell cytoskeleton, we model cytoskeletal attachment and transport of surface molecules by superimposing a biasing factor η on random molecular motion. The probability of successfully diffusing in a particular trial, pdiff, is multiplied by η if the direction of the diffusion move is toward the center, so that η > 1 will favor motion toward the center. It is necessary to show, however, that such an implicit simulation scheme is an accurate simulation of directed transport by the cell’s actin cytoskeleton.
To this end we use a simplified set of simulations consisting of a single N × N lattice, populated by a fixed number M of non-reacting molecules of a single species. Molecules are randomly sampled to undergo diffusion with finite probability pdiff. To simulate attachment to the cytoskeleton, every time a molecule is sampled, a random number trial with probability patt is conducted. If successful, the molecule is considered “attached” to the cytoskeleton. If a molecule that is already attached is selected, it may detach with probability pdet. For an attached molecule, pdiff = 1 if the direction of motion is toward the center, and pdiff = 0 if it is away from the center. Unattached molecules diffuse randomly with probability 0 < pdiff < 1. We run the simulation for a fixed number of time steps and for a range of patt/pdet (α) ratios and compare the distribution of molecules to that obtained from simulations in which the bias factor η is varied and patt = 0 and pdet = 0. As metrics of comparison, we use the mean square displacement (R2) and the mean distance from the center of the lattice (d). Correspondence between values of η and α is established when the values of R2 and d match for the same simulation cutoff time. Shown in Table A1 are the values of R2 and d as η is varied from 1 to 5 and α is varied from 0.01 to 1 for the case of pdiff = 0.1, cutoff time T = 104 time steps, N = 100, and M = 1000 molecules.
From Table A1, we can establish that a bias factor of η = 1.1 approximately corresponds to an attachment/detachment ratio of α = 0.01, η = 1.2 corresponds to α = 0.02, η = 1.5 corresponds to α between 0.05 and 0.1, and η = 2 corresponds to α between 0.2 and 0.5. The above mapping between attachment/detachment ratio and bias factor does not change significantly if we use T = 105 time steps, or change N and M (provided the lattice does not become too crowded).
However, it is important to note that the mapping is dependent on the value of pdiff. For example, for pdiff = 0.5, the values of α listed in Table A1 map to lower values of η (e.g., α = 0.01 corresponds to η = 1.02 instead of η = 1.1). This is due to the fact that as pdiff increases, the probability that a molecule may diffuse away from the center increases more for the case of biased diffusion than for the case of attachment and detachment to the cytoskeleton. For example, for pdiff = 0.1 and η = 1.1, molecules diffuse with probability 0.1 away from the center and 0.11 toward the center, while for pdiff = 0.5 and η = 1.1, molecules diffuse with probability 0.5 away from the center and 0.55 toward the center. By contrast, when η = 1.0 and clustering at the center is achieved by setting patt > 0 and pdet > 0, once a molecule attaches, it will always diffuse toward the center with pdiff = 1 and away from center with pdiff = 0, regardless of the pre-attachment value of pdiff. Thus, while in the case of biased diffusion molecules may diffuse away from the center with higher probability as we increase pdiff, that effect is reduced for the case of attachment and detachment to the cytoskeleton, because the increase in pdiff only affects unattached molecules. As pdiff increases, it becomes “harder” to induce clustering at the center of the lattice by biasing diffusion, and a higher value of η is needed to induce the same degree of clustering. Thus, the same value of α will correspond to a lower value of η as pdiff increases. Similarly, lowering pdiff has the opposite effect, and the values of α in Table A1 map to higher values of η in this case.