Protein Diffusion in Mammalian Cell Cytoplasm

We introduce a new method for mesoscopic modeling of protein diffusion in an entire cell. This method is based on the construction of a three-dimensional digital model cell from confocal microscopy data. The model cell is segmented into the cytoplasm, nucleus, plasma membrane, and nuclear envelope, in which environment protein motion is modeled by fully numerical mesoscopic methods. Finer cellular structures that cannot be resolved with the imaging technique, which significantly affect protein motion, are accounted for in this method by assigning an effective, position-dependent porosity to the cell. This porosity can also be determined by confocal microscopy using the equilibrium distribution of a non-binding fluorescent protein. Distinction can now be made within this method between diffusion in the liquid phase of the cell (cytosol/nucleosol) and the cytoplasm/nucleoplasm. Here we applied the method to analyze fluorescence recovery after photobleach (FRAP) experiments in which the diffusion coefficient of a freely-diffusing model protein was determined for two different cell lines, and to explain the clear difference typically observed between conventional FRAP results and those of fluorescence correlation spectroscopy (FCS). A large difference was found in the FRAP experiments between diffusion in the cytoplasm/nucleoplasm and in the cytosol/nucleosol, for all of which the diffusion coefficients were determined. The cytosol results were found to be in very good agreement with those by FCS.


Introduction
Living cells are multifunctional organisms that exhibit remarkable dynamic phenomena including, e.g., cell motility, and vesicular, cytoplasmic and nuclear transport. The cytoplasm consists of a viscous liquid phase (the cytosol) and a non-liquid phase that will be called here the solid phase. The protein concentration in the cytoplasm has been estimated to be 100 mg/ ml [1], and its total macromolecular concentration (proteins, lipids, nucleic acids, and sugars) can be as high as 400 mg/ml [2]. The cytoplasm can thus be described as a 'molecularly crowded' environment, where macromolecules can occupy 20-30% of its volume [3]. Its solid phase is composed of a dense network of cytoskeletal filaments and membrane structures such as, e.g., the endoplasmic reticulum (ER), Golgi apparatus, and mitochondria [4,5]. Macromolecular diffusion in the cytoplasm can be severely restricted in such an environment [6]. The same applies to the nucleus [7][8][9] that is also composed of a liquid phase, the nucleosol, and a solid phase comprising, e.g., chromatin and proteinaceous nuclear bodies.
Diffusive motion of macromolecules and their binding-dissociation reactions with cellular organelles is a crucial component of cell function, which still need to be clarified. Laser scanning confocal microscopy (LSCM) has become very popular as it allows three-dimensional observation in living cells. LSCM can also be used to perform photo-manipulation experiments such as quantitative fluorescence recovery after photobleaching (FRAP). In FRAP, a region of the cell is exposed to high-intensity laser light, causing the fluorophores within that region to irreversibly lose their ability to fluoresce. Recovery of fluorescence in that region yields information about molecular diffusion and binding in the cell [10,11].
Since the invention of FRAP, several analytical models have been developed to quantify the recovery of fluorescence and thereby diffusion and binding dynamics [12][13][14][15][16][17]. As the internal structure and conditions of the cell are difficult to include in such modeling, several assumptions are made of the system. These assumptions often include infinite, homogeneous fluorophore pools, fast bleaching compared to the time scales of the involved transport processes, and specific shapes of the bleach profiles, conditions that may be difficult to fulfill in FRAP experiments. Models have been suggested that account for diffusion during the bleach phase [15,17], allow for arbitrary bleach profiles [16,18,19], or inhomogeneous distribution of fluorophores inside the cell [20] or of binding sites in the nucleus [21]. Recently the structure of ER [22] has been included when studying protein diffusion in the ER lumen. In these models, however, the constraints imposed by cellular structures have only been partly included.
Fluorescence correlation spectroscopy (FCS) is a method that probes the diffusion coefficient locally, while FRAP probes quite widely the cytoplasm (nucleus). Fluorescence fluctuation microscopy (FFM) combines FCS with LSCM, thereby being able to image the cell environment in which the FCS measurement is conducted. It is expected that the diffusion coefficients measured by FFM and conventional FRAP are quite different unless the internal structure of the cell is included in the latter. There are indeed large differences between the diffusion coefficients reported by FRAP and FFM measurements [14,23,24]. It is evident that, to improve our understanding of dynamic cell functions, better description of macromolecular motion in the cell, and thereby interpretation of experimental results, is called for.
We introduce therefore a completely new approach to protein dynamics in the cell based on describing their diffusive motion in a realistic three-dimensional representation of the cell generated from LSCM data. As it is neither useful nor possible to simulate the dynamics of all proteins of a species in an entire cell at a molecular level of detail, we rely here on mesoscopic methods to model protein distributions instead. This approach is applied here to analysis of FRAP data on the same cell. The model cell takes into account the internal structure of the cell using the inhomogeneous fluorophore distribution at equilibrium. We introduce the new methods by determining the diffusion coefficient of a freely-diffusing model protein, enhanced yellow fluorescent protein (EYFP), in two continuous cell lines, fibroblastlike Norden Laboratory Feline Kidney (NLFK) [25] and cervical carcinoma HeLa [26]. We also determine this diffusion coefficient by FFM for the same cell lines so as to be able to compare the results of the new FRAP analysis with those of FFM analysis. Recent studies using single particle tracking and fluctuation methods have shown that proteins may also undergo anomalous diffusion in the cytoplasm [27,28] as well as in the nucleus [29]. These processes are, however, not considered here as the assumption that all diffusive processes are of Brownian nature suffices to interpret the measured (collective) FRAP data (for theoretical studies of anomalous diffusion see, e.g., [30,31]).

Cell culture
Norden laboratory feline kidney (NLFK) [25] and HeLa [26] cells were grown in Dulbecco's modified Eagle medium (DMEM) supplemented with 10% fetal bovine serum (Gibco, Paisley, UK) at 37uC in the presence of 5% CO2. HeLa cells used for FFM measurements were grown as described in [23]. For live cell microscopy studies, the cells were seeded in 5 cm glass-bottom culture dishes (1.5 thickness, MatTek Cultureware, Ashland, MA). For FFM imaging and measurements, cells were transferred and transfected on 32 mm cover slips as described in [23]. The pEYFP-N3 construct was purchased from Clontech Laboratories Inc. (Mountain View, CA). Transfections were performed with the TransIT-LT1 reagent (Thermo Fisher Scientific Inc, Waltham, MA) according to the manufacturer's protocol. The intracellular localization of the nucleus was visualized by chromatin binding fluorescent histone H2B-ECFP. Cells were transfected with a H2B-ECFP expression vector 24 h after cell seeding.

FRAP experiments
The FRAP experiments were performed on a laser scanning confocal microscope FV1000 with an IX-81 microscope frame (Olympus, Tokyo, Japan) using an Olympus UPLSAPO 606 (NA = 1.2) water immersion objective. The sample stage was heated to 37uC prior the experiments. To image the cell geometry, a confocal stack was acquired before and after the FRAP experiment. The voxel size was adjusted to (200 nm) 3 or (150 nm) 3 . The pinhole size was adjusted to 1 Airy unit. The 514 nm laser line was used for EYFP excitation and the emitted fluorescence was detected using a 530-600 nm band pass filter. Imaging was performed with a laser intensity of 0.1-2 For bleaching a circular (r = 1.85 mm and 2.83 mm) region of interest (ROI) was defined in the middle of the cytoplasm. As bleaching times in FRAP are usually rather large compared to the time scales of the measured diffusion processes, the region of the cell, which is actually bleached, is usually larger than the defined ROI. The size of the actually bleached region and its intensity distribution were measured by bleaching fixed cells (Fig. 1). ImageJ [32] was then used to construct an average shape and intensity profile of that region.
The duration of the bleach process was measured by performing FRAP experiments in which 10 images were collected before the bleach pulse and 1 after the pulse. The bleach time was extracted by measuring the time when the frames immediately before and after the bleach pulse were taken. The average imaging time of one frame was subtracted, and the duration of the bleach process was plotted as a function of iterations (bleaching time), yielding a linear slope. In the LSCM used (Olympus FV 1000), the shortest possible bleach procedure (1 iteration) lasted 36 ms with an additional relay of 18 ms before the next image scan, amounting to 54 ms for the entire process. To achieve enough bleaching for the data analysis, 10 iterations were performed (Fig. 2), and the laser intensity was set to 100% by using an acousto-optical tunable filter.

Image processing
The raw images of the confocal microscope were converted to 8-bit grey scale images. Only linear adjustments of the image brightness and contrast were performed, avoiding saturation. The gray-scale images were colored with an appropriate look-up table and converted to RGB images.

Conventional FRAP analysis
The fluorescence recovery was analyzed in the circular regions described above using the ImageJ [32] and Excel (Microsoft, Redmond, USA) software. Before the measurements, the FRAP data were convoluted with a 363 Gaussian kernel. The data were exported to Excel where their normalization was performed. For normalization two different methods were used. The first normalization (I PM ) used was that of Phair & Misteli [10]: where ROI(t) is the local fluorescence intensity in the bleached region at time t, SROI(tv0)T is the time average of the local fluorescence intensity of the whole bleached region before the bleach pulse, and Cell(t) and SCell(tv0)T are the respective quantities for the entire cell. The second normalization (I A ) used was that by Axelrod et al. [12]: Here ROI(t~0) is the local fluorescence intensity of the bleached region immediately after the bleach phase and ROI(t~?) is the fluorescence intensity after recovery. The time of full recovery can be difficult to determine and requires long imaging times. If the fluorescence is not completely recovered at the end of the experiment, this normalization will lead to too rapid recovery. Therefore we modified the Axelrod normalization (I mA ) such that it could also be used for partial recovery: with ROI(t~end) the fluorescence intensity at the end of the experiment and p the recovery ratio at that time. By definition, the Phair & Misteli normalization always converges to one at full recovery. At partial recovery the value of the Phair & Misteli normalized data can thus be used as the value p needed in the modified Axelrod normalization. The fluorescence intensity data were fitted by the free diffusion model of Soumpasis [13,14]: where t D~r 2 =D f , I(t) is the normalized fluorescence, I 0 and I 1 are modified Bessel functions, r is the radius of the bleached region, and D f the diffusion coefficient of the fluorescent species [33].

FFM measurements
The Fluorescence Fluctuation Microscope (FFM) measurements were conducted with a self made setup in Heidelberg [34]. FFM is a combination of Fluorescence Correlation Spectroscopy (FCS) and Laser Scanning Confocal Microscopy (LSCM). It has an FCS module with a galvanometer scanning unit, attached to the side port of an inverted Olympus IX-70 microscope (Olympus, Hamburg, Germany) equipped with an UplanApo/IR 606 water immersion objective, with a numerical aperture (NA) of 1.2 [34,35].
EYFP was excited with the 488 nm line of an argon-krypton laser from CVI Melles Griot (Bensheim, Germany). The emitted fluorescence from EYFP was recorded between 515 and 545 nm with an avalanche photodiode (APD) (SPCM-AQR-13, PerkinElmer, Wellesley, USA), after passing through appropriate dichroic mirrors and filters for spectral separation and selection. FCS measurements were carried out at laser intensities of 5 to 9 kW cm {2 , and the laser power was adjusted using a polychromatic acousto-optical modulator AOTF Nc (AA Opto Electronic, France). The signals from APD were fed into an ALV-5000/E correlator card (ALV Laser GmbH, Langen, Germany) which recorded the intensity fluctuations and calculated their associated autocorrelation function almost in real time.
The system was carefully calibrated as described in [23] to allow for precise and reproducible measurements.

Construction of the digital model cell and FRAP recovery simulations
For each FRAP experiment we obtained two sets of data: a 3D stack of images depicting the intensity profile of EYFP and H2B-ECFP (histone H2B linked to enhanced cyan fluorescent protein)in the cell before the bleach, and a stack of 2D images depicting a certain cross-section of the cell during the FRAP measurement, 10 frames before the bleach and the rest of the frames from the fluorescence recovery phase. After de-noising the 3D stacks, we used the threshold function of the ImageJ program to segment in the cell the cytoplasm and nucleus using the EYFP and H2B-ECFP stacks, respectively. The nuclear envelope was then generated as a two pixel wide layer between the cytoplasm and nucleus using a self-made code.
The spatial resolution of LSCM, about 200 nm, did not allow segmentation of the more detailed structure. As this 'fine structure' obstructs protein motion and thus affects protein diffusion, we considered the cytoplasm (nucleoplasm) as an 'effective porous medium' that is immobile during a FRAP measurement.
From diffusion in porous media [36] (see Discussion for a simple example), we know that one must distinguish diffusion in the liquid phase (cytosol/nucleosol with D csol =D nsol ) from that in the medium with constrained motion (cytoplasm/nucleoplasm with D cp =D np such that, approximately, D cp~E D csol and D np~E D nsol , where E is the porosity of the medium. As E~E(r), D cp (D np ) is a spatially varying 'effective' diffusion coefficient. A separate diffusion coefficient (D ne ) was assigned to the nuclear envelope described as a permeable membrane. The porosity of the medium was made visible by the heterogeneous equilibrium distribution of fluorophores (proteins), low fluorescence intensity meaning high concentration of 'solids contents', and was thus deduced from the equilibrium fluorescence intensity C 0 (r) (the 3D EYFP stack) such that E(r):C 0 (r)=maxfC 0 (r)g. A 2D cross-section of a typical digital model cell is depicted in Fig. 3.
Due to the low imaging speed of the LSCM used, only a 2D cross section of the cell was imaged during the FRAP experiments. In order to be able to simulate the FRAP recovery in the entire cell, the initial bleach profile had to be extrapolated vertically into the rest of the digital cell. To this end we first determined the relative fluorescence reduction p(x,y) by dividing pixel-by-pixel the first post-bleach image with an average (for noise reduction) of all 10 pre-bleach images of the experimental FRAP data. To enforce the theoretical range of p(x,y) between zero and one, greater valued pixels owing to the noisiness of the experimental data where set to one (flat field correction). The 3D bleach profile was then obtained by multiplying each cross section of C 0 (r) with p(x,y).
The cross section of the cell, which was imaged during the FRAP experiment, was determined by cross-correlating each frame of C 0 (r) with an average of the 10 pre-bleach images of the FRAP stack. The cross-correlation coefficient showed a clear maximum that identified the right cross section.

The lattice-Boltzmann method
The spatial and temporal evolution of diffusion processes is described by the diffusion equation, where r(r,t) is the concentration of diffusing particles and D(r) their possibly locally varying diffusion coefficient. Note that Eq. 5 only accounts for Brownian diffusion processes, in the case of anomalous diffusion, a fractional version of this equation may be used [30]. In Eq. 5 the term is the local diffusive flux of particles. In the present realization we introduce the impermeable solid component in the cytoplasm and nucleoplasm such that an additional flux term, J E , is added to the flux. This term takes care of removal of particles from the nonaccessible regions. The construction of this flux is discussed below. The total flux of diffusing particles is now given by J~J D zJ E , and Eq. 5 can be expressed in the form which is an advection-diffusion equation, where the (local) advective component is given by the additional flux.
In the case of complicated boundary conditions it is, for a numerical realization of Eq. 5, more convenient to start at a somewhat more microscopic level. We thus consider instead the Boltzmann equation [37]. Suitably chosen discrete versions of the Boltzmann equation, in which space, time and velocity are all discrete [37], allow very effective numerical implementations. In the single relaxation time (t) approximation a discrete Boltzmann equation for the distribution function f i (r,t) of particles at point (r,t) moving with velocity v i in the (lattice) direction i, called the lattice-Boltzmann (LB) equation, is given by Here the left-hand side describes the streaming of particles during a time step dt, and the right-hand side models the relaxation of their distribution function towards its local equilibrium, f eq i , on a time scale set by the relaxation time. We have now a three dimensional space and choose a simple cubic lattice with nearest neighbor links only (particles can only move to these nearest neighbors during one time step, which is enough in the case of the diffusion equation [38]). We also allow the particles not to move, and have therefore seven possible velocities (the so-called D3Q7 model [37]) for the particles: i~0, . . . ,6. In this case of an advection-diffusion equation the equilibrium distribution function is given by in which c s is a free numerical parameter (in units of velocity) that determines the proportion of the rest particles, dx is the lattice spacing and w i 's are the D3Q7 weight factors for different discrete velocities: w 0~1 {3c 2 s (dt) 2 (dx) 2 for the rest particles and w i~1 2 c 2 s (dt) 2 (dx) 2 for the other discrete velocities. The second term in Eq. 9 accounts for removing of particles away from the nonaccessible regions. The concentration of particles is given by r(r,t)~P i f i (r,t), and it satisfies (in the continuum limit) Eq. 7 when [38] the diffusion coefficient is given by This diffusion coefficient can be tuned either by changing the relaxation time t, parameter c s , or time step dt. For numerical convenience we fix parameters c s (such that c 2 s~2 7 (dx) 2 (dt) 2 ) and the relaxation time t, and change the diffusion coefficient by tuning the time step.
When applied to modeling a FRAP experiment, the particle density r(r,t) is interpreted as the fluorophore concentration (fluorescence intensity) C(r,t). The additional flux will cancel the diffusive flux into the non-accessible regions filled by membranes, which arises from the concentration gradients in the fluorophore distribution, such that the observed non-homogeneous fluorophore distribution before the bleach will also stay at equilibrium in the simulated model cell. The equilibrium distribution is obtained by setting, at all points r in the cell, Hence the additional flux at equilibrium must be of the form where C 0~C (r,tv0), i.e., it is the (equilibrium) fluorophore distribution before the bleach. During a (FRAP) simulation, the magnitude of this flux at any point in the cell will depend on the actual concentration at that point, which varies in time. Thus, at a given time t, it can be expressed in the form A numerical code was constructed along the lines indicated above, which was capable of simulating the spatial and temporal evolution of the fluorescence intensity in the digital realization of the cell actually measured in the FRAP experiments. When the distribution of fluorescence intensity as measured right after the Figure 4. Visualization of the cross-correlation fitting of corresponding frames. For a given experimental image k, the cross-correlation coefficients c kl (red) each have a global maximum l max (k) (blue crosses). By tuning the parameters t np and t env , l max becomes a linear function of k (black crosses and curve), whose slope determines the simulation time step dt. The deterioration of the maximum in c kl as a function of the experiment frame number stems from the broadening of the bleach profile, which inevitably decreases the relative difference between adjacent frames. Ultimately, this relative difference limits the amount of analyzable experimental frames, which may limit the applicability of the method to slow enough diffusion processes. doi:10.1371/journal.pone.0022962.g004   bleach was taken as the initial condition in the simulation, such simulations could very accurately reproduce the experimentally observed fluorescence recovery.

FRAP data analysis
At every time step during FRAP recovery, the fluorophore distribution was simulated in the whole model cell, and it was recorded in the same cross section as in the measurement. We compared the experimental and simulated frames by cross correlation such that the cross-correlation coefficient was given by Here the subscripts k and l refer to the two series of frames to be compared, v k (x,y) is the pixel intensity of frame k, N is the number of pixels, and v v k and s k are the average intensity and standard deviation of image k, respectively. The experimental and simulated frames could be compared directly using this algorithm. In practice however, the crosscorrelation results were improved greatly if the cell background was removed from all the images (similar to the construction of p(x,y) above), and a mask was used to restrict the analysis to the cell interior. By these manipulations we minimized the perturbing effects of cell motion and deformation.
In the simulations we used three relaxation times to describe the different liquid phases of the cell, the cytosol (t csol ), the nucleosol (t nsol ), and the effective substance of the nuclear envelope (t ne ). Of these three, we fixed t csol for numerical convenience, and the other two were then free parameters. The simulation time step was also a free fitting parameter that eventually determined, together with the values for the three relaxation times, the diffusion coefficients. For a given experimental frame k, the c k,l of Eq. 14 showed a global maximum as a function of the simulation frame number l, denoted by l max (k). The real and digital cells were assumed to correspond to each other when l max (k) was a linear function of k. By varying the values of t nsol and t ne , the linearity of l max (k) was maximized, and the slope of l max (k) directly related the simulation time step dt to the time step Dt used in the experiment (See Fig. 4). These values were then used to determine the values for D csol , D cp (r), D nsol , and D np (r). Finally, we also compared the measured and simulated fluorescence recovery curves for an additional check of consistency.

Validation of methods
In order to test the performance of conventional and the new data analysis methods introduced here, we used the Virtual Cell software [39] to produce data on quasi-2D FRAP experiments with known diffusion coefficients of 10, 25, 40, and 55 mm 2 /s. The bleaching process was modeled as a laser light induced reaction whose creation rate can be described as where L(x,y) is the distribution of laser intensity in the simulation geometry, C EYFP (t) is the concentration of the molecules which are bleached, and V max is the maximum reaction rate. The profile of the laser pulse was either cylindrical with a sharp boundary at a radius of 1.85 mm, or a Gaussian with an HWHM of 1.85 mm (Fig. 5 a). The length of the bleach pulse in the fast and slow bleach simulations was adjusted to 1 ms and 75 ms, respectively. The time lag between the bleach and first recovery data point was 0 ms in the fast and 25 ms in the slow bleach simulations. The simulation time step was set to 0.1 ms or 1 ms, and the fluorescence intensity was recorded at 20 ms intervals. The recovery data were normalized as in the FRAP experiments. We first produced data with different bleach profiles on a geometry (25 mm625 mm, with a thickness of 1 mm and a pixel size of 100 nm) that very nearly conformed to the assumptions made in the Soumpasis method. Next we generated similar data on fluorescence recovery with different bleaching times in two different locations within a 2D digital cell outline determined from an LSCM image of an NLFK cell, with a thickness of 1 mm and a pixel size of 200 nm (Fig. 5 b and c).
As expected, in the ideal case the Soumpasis method recovered the correct value for the diffusion coefficient, especially at the low end of the values used. High diffusion coefficients produced some variation as the size of the periodic box had then a detectable effect on the results. The biggest difference was found for a Gaussian bleach profile placed near the boundary of the cell outline, and a bleach time of 75 ms. In this case the Soumpasismethod diffusion coefficient was 1.56 mm 2 /s, while the correct value was 25 mm 2 /s, a difference by a factor of about 20 (see Fig. 6 for example recovery curves in different experimental conditions).
We then analyzed the same Virtual Cell data by the new LB method introduced here. We found very good agreement in all cases between the correct diffusion coefficient and the one obtained with the LB method, with a maximal deviation of 6%.
We also investigated the effect of data normalization on the FRAP results. In the Soumpasis method data are normalized as in Axelrod et al. [12]. However, the normalization introduced by Phair & Misteli [10], designed so as to include fluorescence loss in the imaging phase, is also very often used. In general, PM normalization increases the diffusion coefficient obtained. This increase seems, however, to be an artefact which arises from the fact that the recovery curve in the Soumpasis method begins at zero and asymptotically approaches one, but when the PM normalization is used, the initial intensity in the recovery phase can be anything between 0 and 1. This problem was demonstrated with the Virtual Cell data for which the Soumpasis method with Axelrod normalization gave the correct result. If we used instead the PM normalization, diffusion coefficients about four times too big were found. A full account of all the analyses done of the various Virtual Cell data is given in Table 1.

FRAP analysis
We performed FRAP experiments on EYFP-expressing NLFK and HeLa cells. When the measured recovery data were analyzed by the Soumpasis method, we found a cytoplasm diffusion coefficient of D~0:75+0:3 mm 2 /s (n = 8) for the NLFK cells and D~1:83+0:28 mm 2 /s (n = 13) for the HeLa cells.
In the new methods introduced, excellent correlations were found between experimental and simulated frames, and the corresponding fluorescence recoveries were also highly consistent (Fig. 7 a). The resulted cytosol diffusion coefficient, D csol , was 55:3+6:8 mm 2 /s for the NLFK (n = 12) cells and 62:2+9:0 mm 2 /s for the HeLa (n = 13) cells. The difference in the liquid phase properties ('viscosity') of these cells is statistically significant (pv0:05), which indicates that they have different macromolecular concentrations. The average cytoplasm diffusion coefficients, SD cp T, were 15:5+2:7 mm 2 /s for the NLFK and 20:6+5:0 mm 2 /s for the HeLa cells, being thus quite similar. In both cell lines the cytoplasm diffusion coefficient, D cp (r), varied significantly (Fig. 7 c).
The emphasis was here on the cytoplasm, but the method automatically produced diffusion coefficients for the nucleus. The nucleosol diffusion coefficients, D nsol , were found to be 28:5+ 16:3 mm 2 /s for the NLFK and 28:2+22:0 mm 2 /s for the HeLa cells, while the average nucleoplasm diffusion coefficients, SD np T, were 18:9+10:8 mm 2 /s (NLFK) and 17:2+13:4 mm 2 /s (HeLa). In these values the uncertainty is obviously rather large as the measurements were not optimized here for their accurate determination. They can, however, be already used for qualitative conclusions.

FFM analysis
As a local measurement technique for gaining further insight into the protein diffusion dynamics, we used FFM (Fig. 8). We measured the cytoplasm diffusion coefficient at 1 to 6 points inside 44 NLFK cells (138 measurement points in total) and at 1 to 9 points inside 50 HeLa cells (198 points in total), from the same cell lines as in the FRAP experiments. The data were fitted by a two diffusing components model [23], whose fast component was estimated to correspond to the cytosol diffusion coefficient determined by the new method introduced here (see the discussion below). In this way we obtained for the cytoplasm an average diffusion coefficient of SD FFM T~60:5+20:2 mm 2 /s for the NLFK and SD FFM T~61:8+19:7 mm 2 /s for the HeLa cells. Note that these values are indeed very similar to the ones found by the new FRAP analysis method for the cytosol, D csol~5 5:3+ 6:8 mm 2 /s for the NLFK and D csol~6 2:2+9:0 mm 2 /s for the HeLa cells.

Discussion
We introduced a new method for modeling protein motion, i.e., evolution of fluorescence intensity in the case of fluorescent proteins, in the entire cell. This method was based on first constructing a three-dimensional digital representation of the cell, which included the cytoplasm, nucleoplasm and nuclear envelope. We furthermore included the effect of internal structures that obstruct protein motion by describing the two cellular compartments as porous media. The equilibrium fluorescence distribution was used to identify the degree by which protein motion is locally obstructed.
The diffusive motion of proteins could then be numerically simulated in a realistic cellular environment. Here we used the lattice-Boltzmann (LB) method for the numerical realization of the diffusion equation. Other methods could also have been used, but the LB method allows an easy implementation of the boundary conditions and an easy generalization to binding-dissociation processes that we intend to include later.
The new modeling instrument was applied to model FRAP experiments that are known to produce typically much lower diffusion coefficients for proteins, when conventional modeling is used to interpret the measured data, than FCS experiments. The fluorescence intensity distribution measured right after the bleach process was used as the starting point for the simulated (postbleach) evolution of that distribution in the entire cell. Thus, in the new method, problems related within conventional modeling to, e.g., intensity normalization, finite volume of fluorophore distribution, and internal membrane structures, were all removed.
As it was expected that the internal membrane structures in the cytoplasm (and nucleoplasm) play an important role in protein transport in the cell, special emphasis was put on properly describing their effect. As described above, such structures could be included as non-accessible regions for protein motion by describing both the cytoplasm and the nucleolasm as porous media. This is not, however, enough. The heterogeneous fluorescence intensity in the cytoplasm/nucleoplasm was interpreted such that it was homogeneous in their liquid phases in which the diffusive motion of proteins only takes place. Distinction was therefore made between diffusive motion in the liquid phases and in the whole cytoplasm/nucleoplasm. The former diffusion coefficients are intrinsic properties of the liquid phases independent of where bleaching is performed. In contrast with this, the cytoplasm/nucleoplasm diffusion coefficient depends on the local membrane structures in and near the region of interest, and varies appreciably.
In order to better understand the distinction between the two types of diffusion coefficient, consider an artificial porous medium a cross section of which is shown in Fig. 9. It can be interpreted as a small region of the cytosol as seen in a confocal microscope image of a cell. The liquid phase is marked blue and the impermeable solid phase is dark brown (its morphology does not represent that of membrane structures in the cytoplasm). We consider diffusion of tracer molecules across the shown structure (homogeneous in the third direction) such that their diffusion coefficient in the liquid phase is set to be 50 mm 2 /s. Diffusion from left to right across the shown medium results in an effective diffusion coefficient of 29 mm 2 /s. This is just a bit smaller that porosity (68%) times the diffusion coefficient in the liquid phase (50 mm 2 /s) because of tortuosity effects (migration paths are in practice longer than the thickness of the region). Inclusion of tortuosity effects is, however, difficult in cellular transport, and they are expected to be rather small on the average. Diffusion of a small non-binding (non-specific binding was assumed to be negligible) fluorescent protein was then analyzed by FRAP, using the new as well as conventional methods, and by FFM, and two different cells were used in these analyzes for generality.
Using FRAP combined with the new analysis method, the cytosol diffusion coefficients were found to be different in the NLFK and HeLa cells, indicating a different macromolecular concentration ('viscosity') in the cytoplasm. In both cases the average cytoplasm diffusion coefficient was about a third of that of cytosol, indicating a similar relative amount of membrane structures in the cytoplasm. FFM analysis of the cytoplasm resulted in diffusion coefficients that were very similar to those found by FRAP for the cytosol. This method probes rather closely the properties of the liquid phase, but if there are 'solid phase' structures near the region analyzed, they however affect [7,8] the result of the measurement. This phenomenon is evidenced by the sizable local variations in the FFM results (in the cytoplasm and in the nucleus [23]). They prevented the detection here of the difference between NLFK and HeLa results. Membrane structures affect the cytoplasm diffusion coefficients in FRAP experiments, and they display strong variation. Evidently it is important to make a distinction in the interpretation of FRAP experiments between diffusion in the cytosol (nucleosol) and in the cytoplasm (nucleoplasm). Using a conventional modeling of the same FRAP experiments, much too low diffusion coefficients were found. The assumptions made in the conventional modeling were not realized in the experimental situation, and no difference was made either between cytosol and cytoplasm diffusion.
Without fine tuning the nucleus and nuclear envelope results we found by the new FRAP method that the D nsol =D csol ratio was about a half in both cells. The nucleosol is thus a more molecularly crowded environment than the cytosol, in agreement with recent results [1]. Both D np =D nsol ratios were about two thirds. It appears that the solid phase affects protein diffusion less in the nucleoplasm than in the cytoplasm [8,23].
For clarity we considered here pure diffusion, but the methods introduced can be extended so as to include interactions of proteins with cellular organelles.