Do membrane proteins cluster without binding between molecules ?

Clustering is a basic event for the initiation of immune cell responses, and simulation analyses of clustering of membrane proteins have been performed. It was claimed that a cluster is formed by the self-assembly induced by protein dimerization with a high binding speed (Woolf and Linderman, Biophys. Chem. 104, 217-227, 2003). We examined the cluster formation with Monte Carlo simulation using two algorithms. The first was that simulation processes were divided into two substeps. All proteins were subjected to movement in the first substep, followed by reaction in the second substep. The second algorithm was that proteins were first selected to react and proteins which did not react were subjected to movement. The self-assembly induced by dimerization was simulated only with the second algorithm. In this algorithm, monomers dissociated from dimers do not move because these monomers are not selected for movement, and a large proportion of such monomers are selected to form dimers in the next step. The self-assembly was again simulated with the first algorithm containing the conditions that monomers dissociated from dimers did not move in the next movement substep. This algorithm seems to be far removed from natural conditions. Thus, it is inferred that the self-assembly induced by dimerization is unlikely in situ, and that some interaction between proteins is required for cluster formation. In contrast to algorithms in previous simulations, our results suggest that it is more appropriate that proteins move to the same direction for a while and reflect when the collision occurs.


INTRODUCTION
Cellular signal transduction is initiated by the binding of a ligand to its receptor.The receptor generally functions in complex forms including homo-and heteromultimers before and after the ligand binding [1][2][3][4][5].Clustering of transmembrane proteins on the cell surface was proposed in the lipid raft model of the plasma membrane [6].Cholesterol, unsaturated sphingolipids and lipid modified proteins etc. do not distribute uniformly in the plasma membrane [7].It is suggested that proteins may exist in "protein islands" connected to the cytoskeleton molecules (protein island model) [8].Foreign antigens are recognized by T cell antigen receptors (TCR) on the cell surface, and the T cells become activated to initiate an immune response [9].The membrane organization of TCR on the T cell surface has been investigated [10][11][12].Similarly, a linker of activated T cells (LAT) was also proposed to exist in "protein islands" on the surface of mast cells and T cells [13].Microscopic techniques have shown separate clusters of TCR and LAT in pre-activated T cells, and these clusters transiently concatenate into microclusters upon antigen recognition [14].The co-stimulation of TCR with CD28 was reported to require co-localization of TCR and CD28 at the plasma membrane [15].
It has remained unclear why such complex formation is required for signal initiation and how the complex is formed.It is hard to answer these questions experimentally because we have still few useful methods to manipulate the complex formation without affecting the function of the proteins themselves.One method to facilitate such examination would be kinetic analysis with the aid of a computer.
Two types of computer simulation techniques are now available: numerical integration of differential equations and Monte Carlo simulation.The former method can address average behavior involving a large number of molecules and stochastic variation.In contrast, the latter can simulate both population behavior and single molecule dynamics.Monte Carlo simulation can also evaluate time-dependent fluctuations involving noise as well as cell-to-cell population heterogeneity [16].
Since one cell contains less than 100,000 molecules of a given membrane protein and because there are variations in biological phenomena, the latter method may be more appropriate.Receptor-ligand formation and clustering of membrane proteins have already been simulated with Monte Carlo techniques [16][17][18][19][20][21], and their results revealed the usefulness of this technique for clarification of biological phenomena.
Various physiological meanings of the clustering of membrane proteins have been proposed [22][23][24][25][26][27], but the mechanism for this cluster formation remains unclear, although a few mechanisms have been proposed [17,18,28].Woolf and Linderman [17] proposed that self-assembly is induced by protein dimerization when the binding speed is higher than the diffusion rate of proteins.We found in this study that different algorithms for Monte Carlo simulation gave different results concerning the cluster formation.The self-organization proposed by Woolf and Linderman [17] was seen in some algorithms, while cluster formation independent of the rate of the dimerization was simulated in other algorithms.We discussed which algorithm was more appropriate for the simulation of complex formation of membrane proteins, and concluded that the self-organization is unlikely in situ.

METHODS
In the present study, a simplified model in which the cell surface is represented as a 2-dimensional plane was assumed, and the cell surface was divided into subspaces.A single subspace was a cubic box with a volume of 166 (5.5 3 ) nm 3 , as described previously [16].One molecule per subspace corresponded to a concentration of 10 mM.Each calculation step was assumed to take 0.02 milliseconds.In all events of our Monte Carlo procedure, real-type pseudo uniform random numbers (N) with the range 0 ≤ N < 1 were generated, as reported previously [29].All proteins were initially distributed into randomly selected subspaces with equal probability.When a selected subspace was occupied, the next subspace was selected randomly.Dimer formation was assumed as follows: The binding of two proteins was accepted when the two proteins occupied neighboring subspaces and N was less than exp(-ΔE1/RT), where ΔE1, R and T are the activation energy, the gas constant and the absolute temperature, respectively.Correspondingly, dimers dissociated when N was less than exp(-ΔE2/RT).
Each protein was assumed to have a movement direction (positive or negative direction on each axis), and a diffusion rate (υ M ).In this study, the movement direction was set randomly, and υ M was set to υ or υ/10, where υ had a Maxwell-Boltzmann distribution from 0 to 999.The probability to have υ (P(υ)) was calculated as fol-

  
When b was set to be 0.005, a Maxwell-Boltzmann distribution of υ was obtained as shown in Figure 1(a), and P(998) and P(999) were 1 and 0, respectively.Proteins moved into their neighboring subspaces according to their movement direction when υ M > τ, where τ is a pseudo uniform random number (0 ≤ τ <1000) obtained by multiplication of N by 1000.Its integer part was used for the rapid simulation.When υ M = 0, the proteins remained in the same subspace.Proteins moved to the opposite side based on periodic boundary conditions when they reached the boundaries of the simulation box.If the opposite side was occupied, the protein was reflected in the mirror direction.If the protein was a part of a dimer, the protein was allowed to pivot around its partner in a random direction.If the target subspace was occupied, the rotation was rejected and not repeated.The movement and rotation of dimers occurred at the same simulation step.
The present simulation included two events; movement and reaction for the formation and dissociation of dimers.We assumed the following methods for the selection of proteins subjected to movement or reaction.MethodA = 0: Simulation processes were divided into two substeps, reaction and movement.All proteins were subjected to reaction and movement in the former and latter substeps, respectively.The reaction substep was carried out after the movement substep.
MethodA = 1: Proteins were first selected to react in each step.Monomers and dimers were converted to dimers and monomers, respectively, according to the reaction probability described above.Proteins that did not react were subjected to movement.
The movement directions of new molecules produced by the formation or the dissociation of dimers were determined randomly, and υ M of such molecules was determined as follows.
MethodB = 0: υ M was always set to 999.MethodB = 1: υ M that had a Maxwell-Boltzmann distribution was set as described above.
MethodB = 2: υ M was always set to 0. When proteins were not reacted, the movement direction and υ M of such proteins were updated as follows.
MethodC = 0: The movement direction and υ M of all molecules were updated immediately before the movement in every step, and molecules moved according to their movement direction and υ M as described above.If the subspace was occupied, the movement was rejected and not repeated.
MethodC = 1: The movement direction and υ M of 0.1 % of molecules selected randomly were updated at every step.If the subspace was occupied, the movement was rejected and not repeated.MethodC = 2: If the subspace was occupied, the molecule was reflected in the mirror direction.υ M was not updated.
MethodC = 3: If the subspace was occupied, the molecule was reflected in the mirror direction, and υ M was updated.
The trajectories of the membrane proteins are shown in Figures 1(b) and (c).
Cluster size was defined as follows: All proteins presented in neighboring subspaces were defined as belonging to the same cluster, and the cluster size was measured by counting all kinds of proteins in the cluster.
The source code of the computing program was implemented using the C-language with Visual Studio C++.net (Microsoft Co.), and the program was run on a personal computer under Windows XP or 2000 (Microsoft Co.).The source code is available from the corresponding author upon request.

RESULTS AND DISCUSSION
In the present simulations, the binding probabilities (exp(-E1/RT)) were set as shown in Table 1.To simulate the binding rate constant (k), the simulation surface was assumed to contain 80 × 80 subspaces and the number of proteins was set to 960.The binding rate constant was calculated from 100 simulated values with 10 different E1 values.The constants were obtained with two other simulation surfaces consisting of 50 × 50 and 100 × 100 subspaces containing 750 and 500 proteins, respectively.The average values are shown in Table 1.The dissociation probabilities (exp(-E2/RT)) were set to one-tenth of the binding probability in all simulations.The average diffusion coefficients calculated from the moving distances of 1000 proteins as described previously [16] are shown in Table 2.In this calculation, proteins are allowed to move even if the target subspace is occupied.In the first simulation (Figure 2), a simulation surface consisting of 80 × 80 subspaces and 960 monomers were set.Fifteen percent of subspaces were occupied by proteins under these conditions.The average cluster size was the same at all binding rate constants in Method 0-1-0 (This means MethodA = 0, MethodB = 1, MethodC = 0) as shown in Figure 2(a).This cluster may be formed without interaction between molecules at a high protein density due to proteins not being distributed uni- Woolf and Linderman [17] proposed that the clustering increased when the binding rate constant was high.In their simulation, molecules were first subjected to reaction and molecules that were not reacted were subjected to movement.The cluster size increased as the binding rate constant increased under their conditions (Method 1-1-0, Figure 2(b)).The same results were obtained in Method 1-0-1, Method 1-1-1, and Method 1-1-5 (data not shown).In this method (MethodA = 1), monomers dissociated from dimers do not move because these The simulation conditions were the same as in Figure 2 except that the cell surface consisted of 50 × 50 subspaces and the number of proteins was initially set to 750.Thirty percent of subspaces were initially occupied by proteins.monomers are not selected for movement.In the next step, a large proportion of such monomers are selected again to form dimers when the binding rate constant is high.This means that a large proportion of dissociated monomers form dimers again without moving at the next step.Therefore, the cluster size is larger at a high binding rate constant.In contrast, when the reaction and movement events are repeated at every step (MethodA = 0), dissociated monomers have the same potential to associate as monomers formed at previous steps.The same results were obtained with different protein densities except that the cluster size and number of dimers increased as the density increased ((a) and (b) in Figures 2-4).To confirm this explanation, the moving energy of dissociated monomers was set to zero in MethodA = 0, i.e., such monomers do not move in the next step.As shown in Figure 2(c), the cluster size increased as the binding rate constant increased.
The increase in the cluster size accompanies a decrease in entropy.In the simple model used in this simulation, no additional energy was supplied for the decrease in entropy when the binding rate constant increased.Therefore, it is reasonable to assume that cluster size is constant at any binding rate constant, suggesting that MethodA = 1 is inadequate.MethodB = 0 and MethodB = 2 seem to be far removed from natural conditions.Therefore, MethodA = 0 and MethodB = 1 seem to be adequate.
The next point is which method is more appropriate in MethodC.In the method described above, the diffusion rates and movement directions of all molecules were updated immediately before the movement in every step.However, it is more reasonable to assume that each molecule has a different molecular activity, namely a different diffusion rate, and keeps the same energy for a while.Therefore, MethodC = 0 is less likely.
The question is thus when does the molecular activity change?We first assumed that 0.1 % of molecules selected randomly were updated in every step (MethodC = 1).The average cluster size decreased as the binding rate constant increased in these conditions, while the decrease in the number of dimers was small (Figure 2(d)).This decrease was similar when the density of proteins increased 2-fold (Figure 3(c)) and small at a low density of proteins (Figure 4(c)).In this simulation, when a protein ran against another protein, its movement was cancelled and its movement direction was not updated.When the binding rate constant was high, each such protein formed a dimer with its neighboring protein immediately, and the movement direction of the dimer was newly assigned.Consequently, the dimer moved away, resulting in a decrease in the cluster size.In contrast, when the binding rate constant was low, proteins that ran against another protein in cluster stayed in the same subspaces for a long time until the proteins were subjected to reaction.This may be the reason for the increase in cluster size at a low binding rate constant.
It is likely that the molecular energy is changed when a collision between molecules occurs in the natural case.In the next simulation, the diffusion direction was updated only when a protein ran against another protein (MethodC = 2).The average cluster size was the same for all binding rate constants (Figure 2(e)).The same results were obtained when both diffusion rate and direction were updated only when a molecule ran into another molecule (MethodC = 3, Figure 2(f)).It is likely that energy is released in a open space even if there is no col-lision.The same results were obtained in MethodC = 4 and 5 which included the updating of 0.1% of molecules at every step (MethodC = 1) in addition to the conditions of MethodC = 2 and 3, respectively (Figures 2(g) and (h)).The cluster sizes were the same again in MethodC = 2 to 5 at the protein densities described in Figures 3 and  4 (data obtained with MethodC = 2 to 4 are not shown).
In MethodC = 1, 4, and 5, 0.1% of proteins selected randomly were updated in every step.It remains unclear whether or not this setting is the most appropriate.The trajectories in MethodC = 5 (Figure 1(c)) were similar to those observed experimentally [34].Although more detailed experimental data are required for more proper setting, 0.1% is probably appropriate.
When the diffusion coefficient of molecules was increased 10 fold (υ M = υ), similar results were obtained except that the cluster size decreased more dramatically as the binding rate constant increased as compared with the lower diffusion coefficient (Figure 5(c)).The diffusion coefficient of membrane proteins observed experimentally was 0.1 to 0.3 μm 2 •sec -1 in prokaryotes [35] and eukaryotes [36,37], and the diffusion coefficient in the setting of υ M = υ was 12.8 μm 2 •sec -1 in MethodC = 1.Therefore, this setting may be less appropriate.It should be noted that the same diffusion coefficient was obtained in MethodC = 1 to 5 because proteins were allowed to move even if the target subspace was occupied when the diffusion coefficient was calculated.
Our present simulations with appropriate algorithms demonstrated that the cluster size was dependent on neither the diffusion coefficient nor the binding speed of proteins at all protein densities tested.Thus, the self-assembly induced by protein dimerization with a high binding speed is unlikely in situ.
GPCRs have been shown to form not only dimers but also oligomers [23,[38][39][40], but structural studies of these receptors have suggested them to have only one protein-protein binding site [41].It may be possible for a membrane protein complex to be formed without such binding site.One possibility is that the hydrophilic surface regions of membrane proteins might bind each other in the membranes.Another possibility is that matrix proteins in the outer or inner cell surface trap membrane proteins in a local area to increase the protein density.It was observed that membrane proteins undergoing Brownian diffusion were confined within a limited area, probably by the binding to a membrane-associated cytoskeleton network [42].In any case, some interaction between proteins may be required for the cluster formation of membrane proteins on the cell surface at a low protein density observed experimentally.

CONCLUSIONS
We examined the cluster formation with Monte Carlo simulation using two algorithms.The first one was that simulation processes were divided into two substeps.All proteins were subjected to movement in the first substep, and then subjected to reaction in the second substep.The second algorithm was that proteins were first selected to react and then proteins that did not react were subjected to movement in each step.The self-assembly induced by protein dimerization with a high binding speed, which was claimed by Woolf and Linderman [17], was simu-lated with the second algorithm, while the cluster size was dependent on neither the diffusion coefficient nor the binding speed of proteins with the first algorithm.In the second algorithm, monomers dissociated from dimers do not move because these monomers are not selected for movement, and a large proportion of such monomers are selected to form dimers before their movement in the next step.The self-organization was again simulated in the former algorithm containing the conditions that the monomers dissociated from dimers did not move in the next movement substep.This algorithm seems to be far removed from natural conditions.Thus, it is inferred that the self-assembly induced by protein dimerization is unlikely in situ, and that some interaction between proteins is required for the cluster formation.
The second algorithm has been used in many previous works, but the present simulation suggests that the first one is more appropriate.We also examined which algorithm was more appropriate for the molecular movement.It has been assumed in many previous simulations that molecules move to the neighboring subspace randomly in each simulation step.In this study, it was shown to be more appropriate that molecules continued to have the same movement direction for a while and the direction was changed at the step selected randomly.Many previous studies adopted the algorithm that the movement was cancelled when the collision occurred.The present study demonstrated that this algorithm was less appropriate, and molecules should change their movement direction in a mirror manner when the neighboring subspace was occupied.It should be clarified in future simulations which interaction is required for clustering of membrane proteins observed experimentally using these appropriate methods.

Figure 1 .
Figure 1.Distribution of diffusion rates and trajectories of membrane protein movement.(a) Distribution of diffusion rates.For details, see text.((b) and (c)) The positions of a given protein were plotted for 1 × 10 5 steps (2 sec) at intervals of 10 steps (0.2 msec).MethodC = 0 (b) and MethodC = 5 (c) were used.The numbers of subspaces and proteins set in this simulation were 80 × 80 and 960, respectively.For details, see text.

Figure 2 .
Figure 2. The average cluster size and the number of dimers when the molecular density was 4950 proteins per μm 2 and υ M was υ/10.The cell surface consisted of 80 × 80 subspaces, and the number of proteins was initially set to 960.Fifteen percent of subspaces were initially occupied with proteins.The diffusion rate of proteins (υ M ) was set to υ/10.Methods used are indicated in the figures.After the reaction reached equilibrium stage, the total number of monomers and dimers was calculated in each cluster at each step, and average values were obtained (closed circles).The number of dimers at each step was calculated and the average percentage of proteins that formed dimers was obtained (open circles).Each point represents the average values obtained from 100 measurements, and standard deviations were less than 5% in all measurements.The horizontal line represents the binding rate constant (k).

Figure 3 .
Figure 3.The average cluster size and the number of dimers when the molecular density was 9900 proteins per μm 2 and υ M was υ/10.The simulation conditions were the same as in Figure2except that the cell surface consisted of 50 × 50 subspaces and the number of proteins was initially set to 750.Thirty percent of subspaces were initially occupied by proteins.

Figure 4 .
Figure 4.The average cluster size and the number of dimers when the molecular density was 1650 proteins per μm 2 and υ M was υ/10.The simulation conditions were the same as in Figure 2 except that the cell surface consisted of 100 × 100 subspaces and the number of proteins was initially set to 500.Five percent of subspaces were initially occupied by proteins.

Figure 5 .
Figure5.The average cluster size and the number of dimers when the molecular density was 4950 proteins per μm 2 and υ M was υ.The simulation conditions were the same as in Figure2except that the diffusion rate of proteins (υ M ) was set to υ.