J. Phys. Chem. B 2008, 112, 143-149
143
Lateral Diffusion of Proteins in the Plasma Membrane: Spatial Tessellation and Percolation Theory Bong June Sung Department of Chemistry, Sogang UniVersity, Seoul 121-742, Republic of Korea
Arun Yethiraj* Theoretical Chemistry Institute and Department of Chemistry, UniVersity of Wisconsin, Madison, Wisconsin 53706 ReceiVed: September 7, 2007; In Final Form: September 26, 2007
The obstructed diffusion of proteins in the plasma membrane is studied using computer simulation and an analysis using spatial tessellation and percolation theory. The membrane is modeled as a two-dimensional space with fixed hard disc obstacles, and the proteins are modeled as hard discs. The simulations show a transition from normal to anomalous diffusion as the area fraction, φm, of obstacles is increased and to confined diffusion for area fractions above the pecolation threshold, which occurs for φm ) 0.22. A Voronoi tessellation procedure is used to map the continuous space system onto an effective lattice model, with the connectivity of bonds determined from a geometric criterion. The estimate of the percolation threshold obtained from this lattice model is in excellent agreement with the simulation results, although the nature of the dynamics in the continuous space model is different from lattice models. At high obstacle area fractions (but below the percolation threshold), the primary mode of transport is a hopping motion between voids, consistent with experiment. The simulations and analysis emphasize the importance of structural correlations between obstacles.
I. Introduction The diffusion of proteins in the plasma membrane is important in many physiological processes, and there is considerable interest in understanding the effect of the structural heterogeneity of the membrane on the dynamics of lipids and proteins.1-3 Extensive theoretical calculations of obstructed diffusion, where the presence of obstacles hinders transport, have been reported.4,5 Although these calculations reproduce the qualitative behavior seen in experiments, there remain quantitative differences. In this work, we study a simple continuous space model for the obstructed diffusion of proteins in the membrane where spatial correlations are present between the obstacles. We show that this model is more appropriate than simple lattice models and can reproduce many experimental observations. We also present an analysis using spatial tessellation and percolation theory. Single molecule experiments provide detailed information on lipid and protein diffusion. These studies have shown that immobile obstacles, for example, junctional complexes bound to the cytoskeleton, can significantly decrease the diffusion constant of proteins in the membrane.6-8 The effect of obstacles is particularly strong in two dimensions (compared to three dimensions) because there are fewer avenues for the diffusing particles to escape from the confinement because of obstacles. The effect of obstacles on protein diffusion is understood at a qualitative level. The presences of obstacles reduces the area available to the diffusion proteins thus decreasing the diffusion constant. As the area fraction of the obstacles is increased, the free area becomes fragmented. There are many paths that the diffusing protein can take which terminate in dead ends and do not contribute to diffusive transport. A consequence is that the available space takes on a fractal nature resulting in anomalous * Corresponding author.
sub-diffusion at intermediate times, where the mean-square displacement, 〈R2(t)〉, is not linear in time, t. Instead 〈R2(t)〉 ∼ tR with R < 1. For long enough times, the diffusion is expected to be normal; that is, 〈R2(t)〉 ∼ t, although this regime may not always be accessible to computer simulations or experiments. For high enough obstacle area fractions, higher than the area fraction at the percolation threshold, the proteins are confined, and 〈R2(t)〉 becomes a constant for large times. The time it takes to cross over from anomalous to normal diffusion increases with increasing obstacle area fraction and becomes infinite at the percolation threshold. Single particle tracking9-11 and fluorescence photobleaching recovery5,12,13 measurements have reported anomalous diffusion of proteins in plasma membrane with R ranging from 0.1 to 0.9. The above picture is particularly clear in lattice models for the plasma membrane. In this model, the obstacles and mobile particles are restricted to the sites of a regular lattice, and particles can move to a neighboring site via a bond provided the target site does not have an obstacle. Theoretical studies14,15 show that the percolation of connected bonds, that is, those that connect two sites neither of which has an obstacle, determines the dynamics. One obtains anomalous, normal, or confined dynamics depending on the fraction, p, of bonds that are connected. Note that p is a measure of the connected bonds and therefore increases when the concentration of obstacles decreases. When p is less than the value at the percolation threshold, there are clusters of connected bonds, but the system is not percolating. In this case, 〈R2(t)〉 ∼ ξ2, where ξ is the correlation length of finite void clusters. At percolation threshold (p ) pc), ξ diverges, and 〈R2(t)〉 ∼ tR at all time scales. For p > pc, anomalous sub-diffusion is observed on short time scales and normal diffusion is observed at long times. The crossover time from anomalous diffusion to normal diffusion decreases as p is increased. A number of computer simulation studies of
10.1021/jp0772068 CCC: $40.75 © 2008 American Chemical Society Published on Web 12/11/2007
144 J. Phys. Chem. B, Vol. 112, No. 1, 2008 lattice models have demonstrated that the presence of obstacles can decrease the diffusion constant and result in anomalous subdiffusion.5,6 Saxton5 reported that R ≈ 0.7 at the percolation threshold and showed that near the percolation threshold there was a broadening of the distribution of diffusion constants and exponents R.5 Saxton has concluded,6 however, that lattice models for obstructed diffusion cannot explain the large decrease in diffusion coefficients seen experimentally, for realistic values of the area fraction of obstacles. Two other classes of theoretical models have been suggested: continuous-time random walk (CTRW)16 and binding models.17 In the CTRW, the solute is assumed to follow a random walk but the probability of taking a step is a powerlaw function of time. The exponent characterizing this waiting time distribution is simply related to the exponent R. This model does not give normal diffusion at long times. Binding models have not received much attention, but in a recent paper, Saxton17 considered the diffusion of particles in a hierarchy of binding sites. With a finite hierarchy, the model gives anomalous diffusion at intermediate times and normal diffusion at long times. The model makes interesting predictions that could be tested via particle tracking experiments. In a series of papers, Kusumi and co-workers2,10 have suggested that a paradigm shift is necessary in our understanding of the transport of proteins in membranes. From high-resolution trajectories in single molecule experiments, they conclude that the plasma membrane is compartmentalized into corrals. The cytoskeleton consists of an array of spectrin filaments (fences) connecting junctional complexes (pickets) and these define the corrals. Within each corral, the proteins undergo Brownian motion and occasionally hop to neighboring corrals. Diffusion is therefore normal at short times with anomalous subdiffusion observed when the mobile proteins sense the presence of the fences and pickets. After a sufficiently large number of hops from corral to corral, the diffusion becomes normal again. Such a behavior is not, however, seen in other experiments.7 In this work, we study off-lattice (continuous space) models of obstructed diffusion. We are interested in addressing several questions raised by recent experiments and theoretical work. Are there significant differences between lattice and off-lattice models? Is it possible to construct an equivalent lattice that reproduces the behavior in continuous space models? Are there simple parameters we can use to characterize the dynamics of proteins in the membrane? Is a strict compartmentalization of space necessary for the hopping motion seen in experiments? We address these issues using a combination of computer simulations and an analysis using spatial tessellation and percolation theory. For this proof of principle study, we employ a crude model for the plasma membrane, that is, a two-dimensional space with obstacles composed of hard discs. Mobile proteins are also modeled as hard discs and can move through void space in a sea of obstacles. The dynamics is investigated using Monte Carlo (MC) simulations, although similar results are obtained using molecular dynamics.18 We use an algorithm based on Voronoi tessellation 18-20 and percolation theory21-23 to map the continuous-space system onto a lattice model, in an effort to characterize the void space and, importantly, the connectivity of void spaces. From the Voronoi analysis, we determine a percolation threshold and find (from the MC simulations) that the transition from confined to normal diffusion occurs at the percolation threshold of the effective lattice. A preliminary communication of some of this work has appeared,18 where MD simulations were used to obtain the dynamic properties.
Sung and Yethiraj There have been previous simulation studies of the diffusion of fluids in random matrices. Chang et al.24 studied the diffusion of hard spheres in a matrix composed of hard spheres (in three dimensions) and showed that the dynamics became heterogeneous as the matrix density was increased. They estimated that a percolation threshold occurred for a matrix volume fraction of about 0.23. Mittal et al.25 presented a simple model, based on estimating the available volume, for estimating the fluid diffusion coefficient. The model was in good agreement with simulations24,25 except near the percolation threshold. There have been theoretical26-29 treatments of the “Swiss cheese“ (or Lorentz) model where the matrix consists of overlapping spheres and the diffusion of a single point particle in this matrix is investigated. Recently, molecular dynamics simulations30 have been used to test the predictions of these theories for scaling exponents near the percolation threshold. The rest of the paper is organized as follows. The model and simulation method are described in section Model and Methods, results are presented and discussed in section Results and Discussion, and conclusions are presented in section Conclusions. II. Model and Methods The model system is composed of fixed (immobile) obstacles and mobile hard discs in two dimensions. The simulation cell is a square of side L with periodic boundary conditions in all directions. The obstacles are hard discs of diameter σm (m for matrix), which is used as the unit of length. The mobile proteins are also modeled as hard discs and have a diameter σf (f for fluid). Unless otherwise noted, σf ) σm ) σ. Initial configurations are generated by sequentially placing obstacles at random locations in the simulation cell. If a trial position for the insertion is overlapped by preexisting obstacles, another random position is tried. This procedure is repeated until a configuration of Nm obstacles is obtained. Then, Nf hard discs of a diameter, σf, are placed in the same way. If there is a percolating free area in the system, particles are inserted into the percolating free area determined by Voronoi analysis. The mobile protein area fraction, φf ( ≡ πσf2Nf/4L2), is fixed at 0.1. The obstacle area fraction, φm ( ≡ πσm2Nm/4L2), ranges from 0.1 to 0.3. In the simulations, 22.4σ eL e126.8σ, Nf and Nm vary between 64 and 6144. For each φm, 10 sets of configurations of obstacles are generated, and the properties reported here are averaged over those realizations of obstacles. To study the effect of the size of obstacles relative to mobile proteins, a range of σf from 0.5 σm to 1.4 σm is investigated for a fixed value of φm ) 0.2. Monte Carlo (MC) simulations are carried out to study the diffusion of mobile proteins. In our previous preliminary communication, discontinuous molecular dynamics (DMD) simulations24 were employed. Molecular dynamics has the advantage that the evolution does not depend on the choice of Monte Carlo moves. The disadvantage, however, is that the proteins move in a ballistic fashion at short times. At long times, however, results from DMD and MC simulations are very similar. In MC simulations, a mobile protein is randomly chosen and translated in a random direction with a maximum displacement, 0.5 σ. If the mobile protein overlaps with other proteins or obstacles, the move is rejected; otherwise, it is accepted. The mean-square displacement, 〈R2(t)〉, is defined by
〈R2(t)〉 ≡ 〈|ri(t) - ri(t ) 0)|〉2
(1)
where ri the position of ith mobile protein, t is time defined as Monte Carlo steps (MCS) and 〈...〉 denotes an ensemble average over both mobile proteins and obstacle configurations.
Lateral Diffusion of Proteins
J. Phys. Chem. B, Vol. 112, No. 1, 2008 145
Figure 1. Schematic of Voronoi tessellation for obstacles (filled circles). Circles are tangent to the three nearest obstacles, and the vertices are the center of the circles. Solid lines are bonds that connect neighboring vertices. bl is the distance between the vertex, l, and the nearest obstacle (radius of the circle plus σm/2), and dm is the distance between the two obstacles that are equidistant from a bond, m.
In Voronoi tessellation, space is tiled into many nonoverlapping convex polyhedra,31 each of which contains one obstacle. A point in space belongs to the polyhedron of an obstacle if and only if the obstacle is the nearest one to that point. To obtain this tiling of space one requires the vertices and bonds that define the polyhedra and these are obtained as follows. For a set of any three obstacles, one can draw a circle that is tangent to all three obstacles. This is done for all sets of three obstacles in the system, and all circles that contain other circles inside them are discarded. This leaves a minimal set of circles, which we refer to as voids. The center of the voids are the vertices of the construction, and bonds connect neighboring vertices (see Figure 1). The bonds and vertices of a polyhedron are therefore equidistant from two and three obstacles, respectively. A vertex is usually intersected by three bonds. The void size is defined as the diameter of the tangential circle, 2bl σm, where bl is the distance between the vertex l and the center of the nearest obstacle. A bond connects two neighbor vertices. The bond is considered as the shortest route of a mobile protein between two voids that the bond connects. The bond width, dm is defined as the distance between centers of two medium particles equidistant from a bond m. After a Voronoi diagram is constructed, the connectivity of each bond is determined using the passage criterion. If a mobile protein of a diameter, σf is too big to pass through the mth bond, that is, dm < σm + σf, the bond is said to be disconnected, and if it is small enough to pass, the bond is considered connected. Once the connectivities of all bonds are found, clusters of connected bonds are obtained using a recursive algorithm. To study the void structure, several properties are calculated, such as the fraction of connected bonds (p), the mean cluster size distribution (S(p,L)) of connected bonds, and the meansquare radius of gyration of a cluster weighted by the number of mobile proteins (〈Rg2〉f). Fraction p is obtained by dividing the number of connected bonds by the total number of bonds. The mean cluster size distribution, S(p,L) is defined as ∞
S(p,L) ≡
sP(s) ∑ s)1
(2)
where s is the number of bonds in a finite cluster and P(s) is the probability that a bond belongs to a cluster with s bonds. Percolating clusters are not included in the definition of S(p,L). 〈Rg2〉f is defined as Nc
〈Rg 〉f ≡ 2
∑ i)1
Rg,i2Nf,i Nf
(3)
where Nc is the number of clusters, R2g,i is the squared radius of
Figure 2. Simulation results for time exponents of 〈R2(t)〉 as a function of φm. The inset shows 〈R2(t)〉 for φm ) 0.1, 0.2, and 0.25 and σf ) σm.
gyration of the ith cluster, and Nf,i is the number of mobile proteins in the cluster. III. Results and Discussion A. Anomalous Diffusion. There are three characteristic regimes in the dynamics of mobile proteins, depending on the area fraction of obstacles. For low values of φm, diffusion is normal, that is, 〈R2(t)〉 ∼ t. As φm is increased, we observe anomalous diffusion, that is, 〈R2(t)〉 ∼ tR at intermediate times with R < 1. In this regime, the diffusion will become normal at long times, but this limit is not accessible in our simulations (or many experiments). For sufficiently high values of φm, above the percolation threshold, the diffusion constant vanishes and the proteins are confined. Figure 2 depicts simulation results for the time exponent R in the anomalous regime, as a function of φm for σf ) σm. Normal diffusion is seen for φm e0.1, the proteins are confined for φm > 0.22, and anomalous diffusion is seen in the intermediate regime. The percolation threshold therefore occurs for φm ≈ 0.22. The inset shows the mean-square displacement for three cases and demonstrates that power-law behavior is seen over several orders of magnitude in time. The exponent R depicted in Figure 2, however, is an effective exponent and reflects the fact that the simulations are for the most part in the crossover regime between anomalous and normal diffusion. The true anomalous exponent is expected to be independent of φm, and R ≈ 0.65 for two dimensions. Our simulation results are in qualitative agreement with previous experiments32 and simulation results.5 Ratto and Longo32 measured the diffusion of proteins in membranes using fluorescence photobleaching recovery experiments and found that R gradually decreased from 1 to 0.7 as the area fraction of the immobile phase was increased. They proposed that the area fraction where R ≈ 0.7 corresponded to the percolation threshold of void space for proteins. B. Void Size Distribution. The diffusion of solutes in matrices is often characterized by the free area available to the mobile species. Even for this simple model of the obstacles, however, the distribution of void sizes is very heterogeneous. Figure 3 depicts voids (open circles) and obstacles (filled circles) for one realization of the matrix with φm ) 0.2. There is a wide distribution of void sizes with some very large voids and some that are too small for the mobile proteins. This distribution is best characterized by plotting the void size distribution, P(σp), which is the probability of finding a void with diameter σp, normalized so that ∫∞0 P(σp) dσp ) 1. Figure 4 depicts P(σp) for several values of φm. As φm is increased, the position of the peak in P(σp) moves to smaller values and the width of the
146 J. Phys. Chem. B, Vol. 112, No. 1, 2008
Sung and Yethiraj
Figure 3. Snapshot of a membrane configuration showing obstacles (filled circles) and voids (empty circles), the latter determined by Voronoi tessellation, for φm ) 0.2.
Figure 4. Void size distribution, P(σp) as a function of the void diameter, σp for several values of φm.
peak decreases. These trends are expected because larger voids become less likely as the area fraction of obstacles is increased. Note, however, that there is no qualitative change in the behavior of P(σp) as the obstacle area fraction is increased above the percolation threshold, while there is a qualitative change in the dynamic behavior; that is, mobile proteins go from being diffusive to being confined at long times. In this regime, the connectivity of the voids is more important than the available free area. C. Connectivity of Voids and Percolation. The Voronoi diagram provides a connectivity map or “road map” for the movement of mobile proteins. Figure 5 depicts a connectivity map for a given configuration of obstacles. In part a, the obstacles are shown with the corresponding Voronoi diagram, and in part b, only the mobile proteins are shown with the same Voronoi diagram. Most of the mobile proteins lie on the solid lines, which implies that the percolating cluster is essentially a road map for the travel of these proteins. The Voronoi diagram presents a means of mapping the configuration of obstacles in continuous space to an effective lattice model. The connectivity properties of this lattice can be characterized by the average fraction of connected bonds, p, and the percolation properties of the lattice for a given value of p. Fraction p decreases as φm is increased, for a given value of σf, because there are fewer paths between voids with width greater than σf. Below some value of p, denoted pc, there is no percolating cluster, and the mobile proteins display confined dynamics; φm,c is the corresponding value of φm. Similarly, p
Figure 5. Connectivity map for φm ) 0.2 and σf ) σ: (a) Snapshot of a configuration of obstacles (red discs) with the corresponding percolating cluster of the Voronoi diagram (solid lines) and (b) mobile proteins in the same configurations (circles of different colors) showing that the Voronoi diagram provides a road map for proteins that are not confined (obstacles are omitted for clarity).
decreases as σf is increased, for a given value of φm, with a percolation threshold occurring for σf ) σf,c. These trends can be seen in Figure 6, which depicts p as a function of φm or σf. The decrease in p is roughly linear with either φm or σf. There is no change in behavior near the percolation threshold, marked with the dashed line p ) pc in the figure. The percolation threshold is obtained via a finite size analysis.22,23 In this method, the mean cluster size distribution S(p,L) is calculated for various values of L. This quantity displays a peak at the percolation threshold, and the peak intensity increases with increasing L although the position is independent of L. As p is increased from zero, the mean cluster size increases and becomes proportional to the system size (diverges in the thermodynamic limit) at the percolation
Lateral Diffusion of Proteins
J. Phys. Chem. B, Vol. 112, No. 1, 2008 147
Figure 6. Average fraction, p, of connected bonds. Open symbols represent the results obtained by changing φm with σf ) σm. Filled symbols represent results obtained by changing σf for φm ) 0.2.
Figure 7. S(p,L) as a function of p for several system sizes. Symbols represent results obtained using the passage criterion with σf ) σm.
threshold. As p is increased beyond pc, however, most of the bonds belong to the percolating cluster, and since percolating clusters are excluded in the definition of S(p,L), S(p,L) decreases. Figure 7 depicts S(p,L) as a function of p for several values of L. The solid symbols represent results using the passage criterion for connectivity discussed so far. A clear peak is seen in S(p,L) which grows as L is increased. From the location of this peak, we determine the percolation threshold to be pc ) 0.526, which corresponds to φm,c ) 0.22 when σf ) σm and σf,c ) 1.1σm when φm ) 0.2. D. Effect of Structural Correlations between Obstacles. There is a significant quantitative difference in the percolation behavior of lattice and continuous space models of obstructed diffusion. Most previous theoretical studies were performed on regular lattices where the obstacles were placed randomly. The two main differences between lattice and continuous-space models are the nature of the lattice, that is, regular versus Voronoi diagram, and the assignment of bond connectivity, which is done using the passage criterion in the continuous space model and randomly assigned in a lattice model. To perform a more quantitative comparison of these models, we investigate a model where the lattice is the same as the Voronoi diagram but the connectivity is assigned randomly using a statistical criterion. For each bond, a random number between 0 and 1 is generated, and the bond is assumed to be connected if the random number is smaller than the desired value of p. This method therefore produces a diagram with the same value of p as in the original Voronoi map but with the bonds assigned in a statistical fashion.
Figure 8. Connected bonds of a Voronoi diagram (φm ) 0.2 and σf ) σ) using (a) passage criterion and (b) statistical criterion. Note that the fractions of connected bonds, p, is the same in the two figures.
The statistical criterion for assigning connectivity results in a percolation threshold at a significantly higher value of p when compared with the passage criterion. At the percolation threshold and with the statistical criterion, the mean number, nc, of connected bonds per vertex is determined solely by dimensionality of space;21 that is, it is independent of the type of lattice. For two dimensions, nc ≈ 2. Since each vertex of the Voronoi diagram is intersected by three bonds, the percolation threshold, which we denote pr, occurs for pr ≈ 2/3.33 Using a method identical to that shown in Figure 7, we calculate pr ) 0.65, which is consistent with the analytical result.33 This is significantly higher than the value of pc ) 0.526. These percolation thresholds are shown as dashed lines in Figure 6. The connectivity of voids is therefore significantly correlated in continuous-space systems even when the obstacles are distributed randomly. This is visually clear when one looks at snapshots of bond configurations for the two cases. Figure 8a depicts a configuration of connected bonds from a Voronoi diagram with φm ) 0.2 and σf ) σ (using the passage criterion).
148 J. Phys. Chem. B, Vol. 112, No. 1, 2008
Figure 9. Mean square displacement, 〈R2(t)〉 for φm ) 0.25 and several configurations of obstacles. Values of 2 〈Rg2〉f are noted in each case. Note that 〈R2(t)〉 is not averaged over media realizations in this figure.
Part b is a snapshot with the same identical configuration of vertices and bonds and value of p but with the connectivity assigned randomly (using the statistical criterion). The lack of correlations in part b results in a less connected diagram which appears to have more bonds, although this is not true. E. Confined Dynamics. For high enough area fractions, there is no percolating cluster, and the mobile proteins display confined dynamics. The mobile proteins are trapped in isolated void clusters, and 〈R2(t)〉 reaches a plateau at long times with the value of 〈R2(t)〉 comparable to 2 〈Rg2〉f, where 〈Rg2〉f is the mean-squared radius of gyration of void clusters weighted by the number of mobile proteins in each cluster. Figure 9 depicts 〈R2(t)〉 (not averaged over realizations of obstacle configurations) for several different configurations. The mean-square displacement shows behavior characteristic of confined dynamics with an initial increase followed by a plateau. Also marked in the figure are values of 2 〈Rg2〉f for the particular configuration of obstacles, and in all cases, the plateau value is comparable to this number. This behavior is consistent with predictions of percolation theory. F. Waiting Time Distribution. An interesting quantity, especially in the context of the CTRW model, is the time a mobile protein spends in a void before leaving the void. We define a waiting time as the time a protein spends in a void (obtained from the Voronoi construction) before it first moves to a different void. The waiting time distribution is calculated by monitoring which void a protein belongs to as a function of time. Note that neighboring voids overlap, and this overlapping region is assumed to belong to the original void the protein was in. Figure 10 depicts the waiting time distribution for several values of φm. The distribution is quite broad, emphasizing the heterogeneous nature of void sizes and connectivity. For example, a protein will take longer to escape from the dangling end of a cluster, compared with a void that is connected to several other voids. At long times, the distribution of a higher φm has a longer tail. The distributions do not display a powerlaw behavior, and the CTRW model is not consistent with these simulations. G. Hopping. For high area fraction of obstacles hopping of the proteins becomes the primary mechanism of transport. In this mechanism, a protein appears to be confined in a void before “hopping“ rapidly to a neighboring void, and the main contribution to translation comes from these hops. Such a form of transport is known to be important in three-dimensional porous materials and glasses, where particles have to overcome an energy barrier to hop to neighboring regions.34 In two-
Sung and Yethiraj
Figure 10. Waiting time distributions for several φm’s, 0.15, 0.2, and 0.25 and σf ) σ.
Figure 11. Coordinate x of a mobile protein as a function of time. Hopping events are indicated by arrows.
dimensional systems, such as the plasma membrane, features out of the plane of the membrane, for example, spectrin filaments attached to the cytoskeleton, are invoked to motivate such motion. We observe clear signs of hopping in this simple continuous space model. Hopping motion is easy to discern from the trajectory of a single particle, although it can be tricky to quantify. Figure 11 depicts the x coordinate of a mobile protein as a function of time for φm ) 0.2. The value of x shows small fluctuations for periods of time punctuated by sharp changes which are the hops. These are marked by arrows in the figure. The trajectory of this protein is plotted in Figure 12. The trajectory has been divided into small time intervals, and a disc is drawn for the position of the protein at each time. Each region where the mobile protein stays before hopping to other regions is drawn in a different color. The trajectory is reminiscent of those seen in experiments. Note, however, that this hopping motion arises solely from the presence of obstacles, and the size of the corrals is determined by the distribution of the sizes of connected voids. IV. Conclusions We study the obstructed diffusion of proteins in a simple model of the plasma membrane using computer simulation and an analysis based on spatial tessellation and percolation theory. The simulations show a transition from normal to anomalous diffusion of proteins as the area fraction of obstacles is increased. The percolation threshold occurs at an obstacle area fraction of 0.22, and for higher area fractions, the dynamics is confined. This percolation threshold occurs at lower area fractions than inferred in lattice simulations. The free area available to mobile proteins is not sufficient to describe the dynamics. In fact, the
Lateral Diffusion of Proteins
J. Phys. Chem. B, Vol. 112, No. 1, 2008 149 These calculations emphasize the importance of continuous space models and structural correlations between obstacles on the dynamics of proteins in the membrane. Clearly, this is a simple model and is not expected to capture all of the phenomena seen in experiments. For example, our model does not consider specific interactions between the proteins and the obstacles or hydrodynamic interactions, both of which could play a significant role. One can envisage extensions of this model where these interactions are taken into account. It might also be possible to dress up the Voronoi lattice with dynamic rules to mimic the hopping of proteins between voids. These extensions would allow for a quantitative comparison between the model and the experiment. Acknowledgment. This material is based upon work supported by the National Science Foundation under Grant CHE0717569. References and Notes
Figure 12. Trajectory of a single mobile protein for φm ) 0.2. Red particles and gray solid lines represent obstacles and a percolating cluster of connected bonds, respectively. Circles represent the position of the protein at different times, and the protein starts from the blue region and ends at the cyan region.
void size distribution functions show little change as one increases the area fraction through the percolation threshold. What does change is the connectivity of the voids. The fragmentation of the free area plays an important role in the dynamics of mobile proteins. A Voronoi tesselation procedure is used to map the continuous space system onto a lattice model composed of voids connected by bonds. A bond between voids is considered connected if it is geometrically possible for a protein to pass directly between the two voids. This procedure provides a road map for diffusion of proteins and is a reliable and accurate method for calculating the percolation threshold. Analysis of the connectivity using the Voronoi map of the continuous space system shows that, even though the obstacles are placed at randomly chosen locations, there are correlations between the connecting paths that make it different from a lattice model with randomly chosen bonds. In fact, the fraction of bonds at the percolation threshold is 0.526 in the continuous space model compared with 2/3 in the lattice model with bonds placed randomly. An interesting feature seen in these simulations is the phenomenon of hopping. At high area fraction of obstacles, the primary mode of transport is a hopping motion between neighboring voids, which is reminiscent of the corral model of Kusumi and co-workers.2,10 The obstacles can therefore be viewed as the membrane proteins that are tethered to the cytoskeleton. Our simulations shed light on some of the questions posed in the introduction. It is possible to construct an equivalent lattice that reproduces some of the qualitative behavior seen in continuous-space models, and the percolation threshold of this model is a useful parameter to characterize the dynamic behavior of proteins. The lattice model is not expected to be faithful to the dynamics, however, as long as all the connected bonds are assumed to be equivalent. A strict compartmentalization of space is not necessary to reproduce the hopping motion seen in experiments.
(1) Glaser, M. Curr. Opin. Struct. Biol. 1993, 3, 475-481. (2) Kusumi, A.; Nakada, C.; Ritchie, K.; Murase, K.; Suzuki, K.; Murakoshi, H.; Kasai, R. S.; Kondo, J.; Fujiwara, T. Annu. ReV. Biophys. Biomol. Struct. 2005, 34, 351-U54. (3) Engelman, D. M. Nature 2005, 438, 578-580. (4) Saxton, M. J. Biophys. J. 1994, 66, 394-401. (5) Saxton, M. J. Biophys. J. 2001, 81, 2226-2240. (6) Saxton, M. J.; Jacobson, K. Annu. ReV. Biophys. Biomol. Struct. 1997, 26, 373-399. (7) Vrljic, M.; Nishimura, S. Y.; Brasselet, S.; Moerner, W. E.; McConnell, H. M. Biophys. J. 2002, 83, 2681-2692. (8) Deverall, M. A.; Gindl, E.; Sinner, E. K.; Besir, H.; Ruehe, J.; Saxton, M. J.; Naumann, C. A. Biophys. J. 2005, 88, 1875-1886. (9) Ghosh, R. N.; Webb, W. W. Biophys. J. 1994, 66, 1301-1318. (10) Ritchie, K.; Lino, R.; Fujiwara, T.; Murase, K.; Kusumi, A. Mol. Membr. Biol. 2003, 20, 13-18. (11) Selle, C.; Ru¨ckerl, F.; Martin, D. S.; Forstner, M. B.; Ka¨s, J. A. Phys. Chem. Chem. Phys. 2004, 6, 5535-5542. (12) Feder, T. J.; Brust-Mascher, I.; Slattery, J. P.; Baird, B.; Webb, W. W. Biophys. J. 1996, 70, 2767-2773. (13) Periasamy, N.; Verkman, A. S. Biophys. J. 1998, 75, 557-567. (14) Havlin, S.; Ben-Avraham, D. AdV. Phys. 2002, 51, 187-292. (15) Ben-Avraham, D.; Havlin, S. Diffusion and Reactions in Fractals and Disorderd Systems; Cambridge University Press: Cambridge, 2000. (16) Montroll, E.; Weiss, G. J. Math. Phys. 1965, 6, 167-181. (17) Saxton, M. J. Biophys. J. 2007, 92, 1178-1191. (18) Sung, B. J.; Yethiraj, A. Phys. ReV. Lett. 2006, 96, 228103. (19) Kerstein, Phys, A. R. J. A: Math. Gen. 1983, 16, 3071-3075. (20) Elam, W. T.; Kerstein, A. R.; Rehr, J. J. Phys. ReV. Lett. 1984, 52, 1516-1519. (21) Kirkpatrick, S. ReV. Mod. Phys. 1973, 45, 574-588. (22) Stauffer, D. Phys. Rep. 1979, 54, 1-74. (23) Isichenko, M. B. ReV. Mod. Phys. 1992, 64, 961-1043. (24) Chang, R.; Jagannathan, K.; Yethiraj, A. Phys. ReV. E 2004, 69, 051101. (25) Mittal, J.; Errington, J. R.; Truskett, T. M. Phys. ReV. E 2006, 74, 040102. (26) Halperin, B. I.; Feng, S.; Sen, P. N. Phys. ReV. Lett. 1985, 54, 2391. (27) Feng, S.; Halperin, B. I.; Sen, P. N. Phys. ReV. B 1987, 35, 197. (28) Rintoul, M. D.; Torquato, S. Phys. ReV. E 1995, 52, 2635. (29) Machta, J.; Moore, S. M. Phys. ReV. A 1985, 32, 3164. (30) Ho¨fling, F.; Franosch, T.; Frey, E. Phys. ReV. Lett. 2006, 96, 165901. (31) Okabe, A.; Boots, B.; Sugihara, K.; Chiu, S. N. Spatial tessellations; John Wiley and Sons: New York, 2000. (32) Ratto, T. V.; Longo, M. L. Langmuir 2003, 19, 1788-1793. (33) Hsu, H.; Huang, M. C. Phys. ReV. E 1999, 60, 6361-6370. (34) Jagannathan, K.; Sung, B. J.; Yethiraj, A. Phys. ReV. Lett. 2006, 97, 145503.
