Three-dimensional dynamic deformation of a red blood cell in a dual-trap optical tweezers is computed with the elastic membrane theory and is compared with the experimental results. When a soft particle is trapped by a laser beam, the particle is deformed depending on the radiation stress distribution whereas the stress distribution on the particle in turn depends on the deformation of its morphological shape. We compute the stress re-distribution on the deformed cell and its subsequent deformations recursively until a final equilibrium state solution is achieved. The experiment is done with the red blood cells in suspension swollen to spherical shape. The cell membrane elasticity coefficient is obtained by fitting the theoretical prediction with the experimental data. This approach allows us to evaluate up to 20% deformation of cell’s shape
© 2010 OSA
Optical tweezers has been used extensively for manipulating living biological cells [1–3]. This approach allows us to impose systematically controlled stress states on the cells to study the constitutive response of the cell membrane and cytoskeleton to physical, chemical, and biological environmental conditions. Various multiple-trap optical tweezers have been proposed for manipulating the cell and measuring the elastic properties of an isolated cell. The optical radiation forces can be applied via the beads attached to the cell [4,5]. In this case the mechanical loading is localized in the regions where the beads are bound to the cell. The multiple traps have been also applied directly to the cells [6–8] without any physical contact. For example, Bronkhorst et al used three-traps to bend the discotic red blood cells (RBCs). In this experiment, the three trapping focal spots were first aligned in the interior of the biconcave shape of a RBC and then the middle focal spot was moved normally to the line connecting the other two spots to deform the cell into a parachute shape. As the experiment was to measure the cell’s recovery time, there was no need to evaluate the radiation stress and the corresponding deformation of the cell in . Liao et al stretched the bi-concave RBC directly by a trapping beam, which was jumped over a distance in the interior of the cell. The cell’s elongation was measured as a function of the beam jumping distance. The experimental result was explained qualitatively by the 2D stress distribution on RBC discotic disk .
In this paper, we investigate the deformation of a cell in the dual-trap optical tweezers by a quantitative analysis of the 3D radiation stress distribution and the corresponding 3D deformation for a RBC swollen to spherical shape. Deformation of a spherical RBC under radiation pressure has been studied experimentally and theoretically with a two-counter-propagating-beam laser stretcher [9,10]. In the case of the dual-trap tweezers , however, the two parallel trapping beams focused in the interior of the cell are laterally separated by a distance. This destroys the axial symmetry and therefore increases the complexity of the computation. Moreover, the previous works considered only the static deformation [9,10], that is to compute first the radiation stress on a spherical cell and then the consequent cell’s deformation using linear elastic theory of membrane. This analysis is valid only for small deformation of about 5-10%. In fact, the RBCs usually undergo much larger deformation. When the cell is deformed under initially applied stress, the radiation stress will be redistributed on the deformed shape if the trapping laser beams are not removed. Then, the cell will be further deformed under the redistributed stress. In this paper we consider this dynamic regime of the deformation and compute the stress redistributions and the subsequent deformations in a recursive manner using the finite element method, until a finale equilibrium state is reached. Fitting the theoretical prediction with the dynamic deformation regime to our experimental data shows that the proposed approach can be valid for cell’s deformation up to more than 20%. The fitting gives the elasticity of the human RBC membrane.
2. Stress distribution
Assume a RBC, swollen initially to a spherical shape, is centered in the Cartesian coordinate system. Two trapping beams coming along the + z axis are shifted along the ± x axis and focused at (x = ± D/2, y = z = 0), respectively. As the RBC diameter 2ρ ~7 µm is larger than the wavelength λ = 1.06 µm, the ray optics approach was used to compute the radiation pressure approximately [1,3,9]. According to the law of momentum conservation, the radiation stress isEq. (1). As n2/n1>1, the stress is negative with respect to the incident rays, resulting in a stretching force pointing outward from the cell. The cell will be trapped to an equilibrium position where the net force on the cell is equal to zero. In this position, the cell’s center will not remain in the plane of the focus z = 0. However, we ignored this slight offset in the calculation of the stress and the deformation, for the sake of simplicity. The mathematical formulas for calculating the stress on the sphere by a shifted focusing beam using the ray-tracing approach are given in Appendix A.
In an alternative approach we used the fact that the local radiation force is normal to the interface [3,9], as the optical reflection and refraction do not change the tangential component of the photon’s momentum to the interface, so that so that we considered only the incident photon momentum component normal to the cell surface. In this case, , , and were all aligned to the normal of the surface. The calculation still follows what described in Appendix A, but the components Qx, Qy and Qz in Eq. (A1) were not required. The Fresnel reflection coefficients are still computed with Eq. (A2). The two approaches gave identical result for the spherical cell. However, the latter approach is simpler. Moreover, when the cell is deformed, the vector in Fig. 6 will be no longer normal to the surface and Eq. (A1) is no longer valid, so that the latter approach must be used with the replaced by the normal of the deformed surface for computing the incident and refraction angles used in Eq. (A2).
The stress on the rear surface of the sphere was computed by taking the transmitted ray as the incident ray with the intensity reduced by T and neglecting the absorption by the cytosol. As n2-n1 ~0.043 is very small, the reflectance R at the cell/buffer interface is small, we neglected the radiation stress associated with the third and subsequent reflections inside the cell.
We computed the radiation stress on the spherical RBC surface associated with the two trapping beams focused at (x = ± D/2, y = z = 0), respectively and then added the two stress distributions. The results are shown in Fig. 1 When the beam separation D = 0, the stress distribution is revolutionary symmetric with respect to the z-axis, and there is a crown-shaped region around the equator of the sphere in the x-y plane, on which only few or no optical rays incident, so that the stress is low or equal to zero, as shown in Fig. 1(a). The width of this region depends on the numerical aperture of the trapping beam and the size of the cell. We then shift the beam focus along the x-axis. We observed that the stress distribution remains practically unchanged in the range of D < 2 μm. In the view of projection along the z axis the elongation of the cell along ± x axis is negligible for D < 2 μm. This is in agreement with the experimental observation  and with the photographs shown in Fig. 5(a) and 5(b) for D = 1,27 – 2,54 μm. When D continues to increase the stress is concentrated gradually towards the x-y plane and towards the ± x directions with the magnitude of the stress peaks increasing significantly, as shown in Fig. 1(b)–1(d).
The trapping field was also solved by the generalized Lorentz-Mie scattering theory (GLMT) with the T-matrix  and by the numerical FDTD algorithm (OptiwaceTM). Then, the radiation stress was computed using the Maxwell stress tensor. As shown in Fig. 1(d)–1(f), the ray-tracing, the T-matrix and the FDTD method produced similar results with the computed peak stress at the first surface of the cell as 0.979, 0.923 and 0.980 N/m2, respectively, for the beam separation D = 5.07 μm. The differences in the stress values are mostly due to the beam intensity distributions and the beam models. In fact, in the ray tracing, the trapping beam was modeled as a bundle of the focused optical rays with the Gaussian intensity distribution. In the FDTD a paraxial Gaussian beam model was used. On the other hand, in the GLMT the trapping beam was modeled as the TEM00 (or LG00) beam satisfying the vector Helmholtz equation for the strongly focus beam and was then expanded by the multipole expansion into the orthogonal set of vector spherical wave-functions (VCWF) taking into account the high NA. Moreover, the expansion coefficients were determined by matching the far field distribution of the VSWF expansion to the far field Gaussian beam .
The calculations with both the ray-tracing and T-matrix showed that the stress peak values on the rear surface (z > 0) is about 10% higher than that at the front surface (z < 0). This is because Eq. (1) is based on Minkowski model according to which the photon momentum increases with the refractive index n . For instance, we evaluate the stresses at the front and rear surfaces with n1 = 1.335 in the buffer and n2 = 1.378 inside the cell and the reflectance and transmittance, R and T, at the normal incidence as:
As n2> n1 and R > 0, we have σrear > σfront. However, the difference in σrear and σfront has no direct consequence on our experiment, as only the cell’s elongation along the ± x axis will be taken as a measure for cell’s deformation. In Section 4 we will choose the ray tracing approach and embed its MatlabTM codes into the Comsol multiphysicsTM modules for computing the dynamic deformation with the finite element method, for the sake of simplicity.
3. Static deformations
In the well-known liquid-interior model, the RBC consists of thin membrane containing incompressible fluid (cytosol) [14,15]. The human RBC membrane consisting of a phospholipid bilayer, a spectrin filament network and transmembrane proteins can be treated as a hyperelastic effective continuum layer. The ratio between the thickness h of the membrane and the cell radius ρ, h/ρ ~1% . The membrane has finite flexural stiffness to resist in-plane tensile or compressive forces, but cannot support off-plane bending and twisting forces. In our experiments the RBC samples were swollen from physiological biconcave shape to a spherical shape by the osmotic pressure. Once captured by the dual trapping beams, the RBC received radiation pressure as an additional external load and was deformed to a new state of equilibrium in which the external load is balanced by the membrane internal force. For a spherical membrane of radius ρ with only the radial stress σr(φ,θ), the equilibrium equations are given by .Eqs. (2a)–(2c), and then associate the stresses to the local strains prescribed by Hooke’s law in the limit of the linear elasticity. The deformation of the cell is described as displacement vectors of the material points on the membrane and can be solved from the local strain by the constitutive equations, which are for the spherical membrane
In the case of spherical cell the equilibrium equations and the constitutive equations take simple forms as Eqs. (2) and (3). However, in the dual-beam optical tweezers the external load σr = σr(φ,θ) is not symmetric to the z-axis due to the lateral shift of the two beams. In this case we can still solve the equilibrium and the constitutive equations analytically by using the Fourier expansions with respect to θ. The mathematical detail is given in Appendix B.
4. Dynamic deformation
The static solution is valid only when cell’s deformation is small. When the beam separation D is large and the cell’s deformation is significant, the re-distribution of the radiation stress on the deformed surface needs to be considered. Furthermore, the stress re-distribution will lead to cell’s re-deformation, which in its turn will lead to new re-distribution of the stress. The process continues until a final equilibrium state is reached for a given incoming optical field. In this case the equilibrium and the constitutive equations for the deformed cell can no longer take the simple forms as Eqs. (2) and (3) and analytical solutions would be more difficult. We used the finite element solution in ComsolTM structured mechanics module to compute the linear elastic deformations of the deformed cell. To implement the computation we embedded our MatlabTM codes for computing the radiation stress with the ray optics approach to the ComsolTM module.
We then input the spherical RBC of radius ρ, the membrane thickness h, a test value of the elasticity coefficient E and the beam separation value D, and we selected a uniform deformable moving mesh on the cell surface. As in our model the net radiation force applied to the cell is not zero, the cell’s position would be shifted in the computing. However, this shift is not of interest when computing cell’s deformation. We set the constraints that the 6 poles of intersection of the cell surface with the x, y and z axes can move only along the x, y and z axis respectively for preventing the cell from shift in the space. We launched the stationary solution solver in the ComsolTM. The stress distribution was first computed by our MatlabTM codes, and then the cell’s deformation under the computed stress distribution was computed by a linear solver SPOOLES. Once SPOOLES completed its execution, the cell was deformed along with the moving mesh, so that the stresses on the deformed mesh nodes were re-computed using the embedded MatlabTM codes. Then, the SPOOLES was launched again to compute the deformation under the new stress distribution of the deformed mesh nodes. The iterative process continued until a convergence of the solution was finally reached after 4-5 iterations, and the deformed cell reached a final equilibrium state. We found that the cell’s deformation computed in the first iteration of the linear solver SPOOLES with the stress distribution on the initially spherical cell is identical to the static deformation computed by the analytical method described in the Appendix B. This agreement validates our ComsolTM solution.
Figure 2 shows the stress redistribution on the deformed cell. With the beam separation D = 2.54 μm and the radius ρ = 3.86 μm the stress distribution on the initially spherical cell is shown in the top-left in Fig. 2, under which the spherical cell was deformed and the stress redistribution on the deformed cell was computed in the iterations as shown in the clock-wise order in Fig. 2. We see that the peaks of the redistributed stress became wider and the peak values decreased as the cell was deformed progressively. As a result, the final cell’s deformation in the dynamic regime is quite different from that computed in the static regime without considering the re-distribution of the stress. The dynamic regime is closer to the physical reality for larger deformations.
Figure 3 shows the computed 3D deformation of a spherical RBC in the dual beam optical tweezers as a function of the beam separation distance D, which are now ready to fit the experimental results.
A schematic diagram of the dual-trap optical tweezers setup is illustrated in Fig. 4 . A linearly polarized laser beam (λ = 1064nm, 1W, Nd: YVO4 cw laser) for optical trapping passed through a beam expander (BE) which expanded and collimated the beam such that the beam diameter overfilled the back aperture of an oil immersion microscope objective (NA = 1.25, 100X). The laser beam was split into two via a polarizing beam splitter (PBS) and recombined by a second polarizing beam splitter (PBS) to form two parallel trapping beams, focused in the focal plane of the objective lens; the optical powers of the two beams can be adjusted byrotating the half-wave (HW) plate. One of the laser beams passed through a telescope consisting of two lenses with 200 mm focal length. The first lens was mounted on a motorized translational stage (MS) for displacement in the direction transversal to the optical axis. The telescopic imaging arrangement transformed the transverse displacement of the first lens into a lateral shift of the focal spot of the trapping beam without any beam walk-off at the entrance pupil of the objective lens. The RBC was trapped and stretched in the focal plane. Wide-field images of the RBC were captured by a CCD camera for observation and analysis.
The RBC sample was swollen to spherical and suspended in saline (150mOsm) buffer solution. The distance D between the focal spots of the two trapping beams was increased from 0 μm to 6.34 μm in 6 discrete steps. As the cell’s reaction to the external load is typically slow, in the experiments we increased the beam separation D stepwise and always allowed enough time for the cell to reach equilibrium before the next increase of D. In each step we measured the length of the major axis of the RBC along the direction of the beam separation with a pre-calibrated length-scale on the CCD image. The elongation was measured for the cell without drug treatment and with 1mM N-ethylmaleimide (NEM) treatment for 30 minutes as a function of the dual beam separation distance D; Experimental micrographic image examples are shown in Fig. 5(a) and 5(b) for the RBCs without drug treatment and with 1mM NEM treatment for 30 minutes, and for the different dual beam separation distances respectively.
From the theoretical predictions on 3D cell’s deformation shown in Fig. 3, we see that the projections of the 3D shapes of the deformed cells along the z-axis to the x-y plane are the ellipses elongated in the ± x directions as observed in the experiments, and shown in Fig. 5(a) and 5(b). The cell’s deformation along the z-axis was not readily observable in the experiments. The fitting parameters were the product of the membrane elasticity and the membrane thickness Eh and the cell’s initial radius ρ. As ρ can vary in the experiments among the group of 30 samples, we took ρ as a fitting parameter. We found Eh = 15.99 (μN/m) and ρ = 3.554µm for the normal RBC’s and Eh = 24.39 (μN/m) and ρ = 3.697µm for the RBCs with the drug treatment. The corresponding shear modulus are Gh = Eh/2(1 + v) = 5.33 and 8.13 (μN/m), respectively.
Our fitting is valid for more than 20% of deformation, measured as a ratio of the elongation of the long axis length to the cell’s initial diameter. A RBC in suspension can usually restore its initial shape after such a large deformation . This is similar to a balloon in air, which can undergo a large deformation keeping in a constant volume before explosion. Note that the shape deformation can be caused by some simple displacements and some deformations of the membrane parts. When the RBC is elongated by 20% the deformation (stretching or compression) of the membrane parts can be much less than 20% and therefore still in the linear elastic range.
We have implemented theoretical calculation on the dynamic deformation of the red blood cell with the finite element method and Comsol MultiphycisTM modules along with embedded codes to compute the stress distribution on the deformed cell. Our tool thus combines the geometrical optics and the structure mechanics into the single calculation. This allows computing the stress re-distribution and cell’s re-deformation in the iterations. Our theoretical predictions fit to the experimental data permitting differentiating the red blood cells with and without the drug treatment and obtaining the cell membrane’s elasticity coefficients. This approach allows us to evaluate large deformation of the cells in the dual-trap optical tweezers.
Consider one of the rays in the focused trapping beam incident at an arbitrary point on the sphere B(r, φ, θ) and refracted by the cell/buffer interface toward the focal point F (x = D/2, y = z = 0), as shown in Fig. 6. Take the refracted ray, which spans the incident plane with the normal to the surface at point B. In this plane the refracted angle may be computed, according to the Snell’s law, as and the incident angle is δ = sin−1 (n2sinβ/n1). Then the vector can be written as
where the angle ϕ between and the x-axis is given by ϕ = tan−1((z2+y2)1/2/x). The Fresnel reflection coefficients are computed for the polarization normal to and parallel to in the incident plane, respectively, as
As the polarization of the trapping beams was not controlled in the experiments, the mean intensity reflectivity and transmittance T=1-R were used in Eq. (1). The nominal power P in Eq. (1) of the incident ray in the Gaussian beam is determined in the parallel plane to the x-y plane that sections the trapping beam at the incident point B, , where the Gaussian beam width w depends on the distance z from the plane of section of the beam to the focal spot. Thus the stress σ(φ,θ) can be computed from Eqs. (1) and (A1).
That results in one differential equation for the Fourier coefficient Nm as
which is solved by Bernoulli method as
where Am and Bm are integration constants. Sum and subtract Nm and Lm allow us to obtain
where Cm and Dm are integration constants
Finally, from Eq. (2c) we have
The deformations (u, υ, w) are reconstructed by the Fourier series with respect to θ using (um, υm, wm). The integration constants Am, Bm, Cm, and Dm are useful to reconstruct the function A, B, C and D, which were trimmed in order to keep the cell’s volume constant during the deformation as the liquid inside is incompressible.
References and links
3. J. Guck, R. Ananthakrishnan, H. Mahmood, T. J. Moon, C. C. Cunningham, and J. Käs, “The optical stretcher: a novel laser tool to micromanipulate cells,” Biophys. J. 81(2), 767–784 (2001). [CrossRef] [PubMed]
4. M. Dao, C. T. Lim, and S. Suresh, “Mechanics of the human red blood cell deformed by optical tweezers,” J. Mech. Phys. Solids 51(11-12), 2259–2280 (2003). [CrossRef]
5. S. Hénon, G. Lenormand, A. Richert, and F. Gallet, “A new determination of the shear modulus of the human erythrocyte membrane using optical tweezers,” Biophys. J. 76(2), 1145–1151 (1999). [CrossRef] [PubMed]
6. P. J. H. Bronkhorst, G. J. Streekstra, J. Grimbergen, E. J. Nijhof, J. J. Sixma, and G. J. Brakenhoff, “A new method to study shape recovery of red blood cells using multiple optical trapping,” Biophys. J. 69(5), 1666–1673 (1995). [CrossRef] [PubMed]
7. W. G. Lee, H. Bang, J. Park, S. Chung, K. Cho, C. Chung, D.-C. Han, and J. K. Chang, “Combined microchannel-type erythrocyte deformability test with optical tweezers,” Proc. of SPIE.6088, 608813–1-12, (2006).
8. G. B. Liao, P. B. Bareil, Y. Sheng, and A. Chiou, “One-dimensional jumping optical tweezers for optical stretching of bi-concave human red blood cells,” Opt. Express 16(3), 1996–2004 (2008). [CrossRef] [PubMed]
10. J. Guck, R. Ananthakrishnan, H. Mahmood, T. J. Moon, C. C. Cunningham, and J. Käs, “The optical stretcher: a novel laser tool to micromanipulate cells,” Biophys. J. 81(2), 767–784 (2001). [CrossRef] [PubMed]
11. T. A. Nieminen, V. L. Y. Loke, A. B. Stilgoe, G. Knöner, A. M. Brańczyk, N. R. Heckenberg, and H. Rubinsztein-Dunlop, “Optical tweezers computational toolbox,” J. Opt. A, Pure Appl. Opt. 9(8), S196–S203 (2007). [CrossRef]
12. T. A. Nieminen, H. Rubinsztein-Dunlop, and N. R. Heckenberg, “Multipole expansion of strongly focussed laser beams,” J. Quant. Spectrosc. Radiat. Transf. 79–80, 1005–1017 (2003). [CrossRef]
14. F. Xu, J. A. Lock, G. Gouesbet, and C. Tropea, “Optical stress on the surface of a particle: homogeneous sphere,” Phys. Rev. A 79(5), 053808 (2009). [CrossRef]
15. Y. C. Fung, “Theoretical considerations of the elasticity of red cells and small blood vessels,” Fed. Proc. 25(6), 1761–1772 (1966). [PubMed]
16. P. B. Bareil, Y. Sheng, Y. Q. Chen, and A. Chiou, “Calculation of spherical red blood cell deformation in a dual-beam optical stretcher,” Opt. Express 15(24), 16029–16034 (2007). [CrossRef] [PubMed]
17. E. Ventsel, and T. Krauthammer, Thin plate and shells (Marcel Dekket, New York, 2001).