Three-dimensional light-scattering and deformation of individual biconcave human blood cells in optical tweezers

Lingyao Yu; Yunlong Sheng; Arthur Chiou

doi:10.1364/OE.21.012174

1. Introduction

Optical tweezers are useful to study deformation characteristics of a single living human red blood cell (RBC). One way to manipulate an RBC in the optical tweezers is attaching two silica microbeads to the cell and using two laser beams to trap the microbeads and stretch the cell [1]. This method has the advantage of imposing a simple stress state, such as the tensile stretching, to the cell with controllable amount of deformation. A potential disadvantage is that the positions and areas of the contact regions between the cell membrane and the microbeads may be difficult to control with high accuracy in the experiments, resulting in uncertainty in the stretching force calibration. Direct trapping and manipulating a cell by the optical tweezers without mechanical contact to the cell was first demonstrated by Ashkin in 1987 [2]. However, when a naturally biconcave-shaped RBC is trapped by a linearly polarized focused beam, the cell will be twisted and rotated in the isotonic solution [3]. Therefore, two parallel trapping beams were used to stably trap a biconcave RBC with controllable cell orientation [4–6]. In this case the trapped biconcave RBC will be flipped to a stable orientation with its platelet parallel to the plane defined by the two trapping beam axes [4]. As a consequence, the RBC is always oriented such that the trapping beams enter from the side of the RBC with the laser beam axes parallel to the biconcave platelet. More recently, the dual-trap optical tweezers were used in the micro-rheology applications with the two beams trapping the RBC at its opposite edges, with or without attaching microbeads. With one of the beams oscillating at a small amplitude and a frequency up to 1 KHz to introduce fluctuation, the phase signal of the cell edge displacement vibration was recorded as a function of the oscillation frequency to probe the RBC viscoelastic properties [7, 8].

Modeling and numerical simulation is important for understanding the physical phenomena and interpreting the experimental data. When the trapped particle is a deformable cell, such as the RBC, the stress distribution on the cell is first computed to evaluate the consequent cell’s deformation. The stress distribution on a biological cell is usually computed by the geometrical optics ray tracing and the law of photon momentum conservation. The complexity in the stress calculation increases with the complexity of the cell’s morphological shape [9, 10]. The geometric optics ray tracing is more difficult to perform for the biconcave-shaped RBC. It is even more difficult to simulate the final equilibrium state of the trapped cell. One has to compute the final stress distribution on the deformed cell. In this case, we prefer to compute the scattered field by the trapped cell and to obtain the stress distribution via the Maxwell stress tensor. In fact, the scattering of light by the RBC has been an important research topic with applications to disease diagnostics and drug therapies [11–14]. For example, the scattering technique has been used for measuring the morphological parameters of clinical relevance of the cell by inverse scattering methods [15, 16]. Light-scattering by the RBC has been computed with the discrete sources method (DSM), discrete dipole approximation (DDA), T-matrix, finite difference in time domain (FDTD) and many other methods [17]. However, most works deal with neither side-incidence nor focused trapping beams [18, 19].

In this paper, we report the full 3D simulations and modeling with the finite element method (FEM) for the scattering and the deformation of the biconcave RBC in the dual-trap optical tweezers. We computed the scattering of a trapped biconcave RBC with the FEM electromagnetic module, the stress distribution via the Maxwell stress tensor and the cell’s deformation with the FEM mechanics modules. Furthermore we iteratively coupled the two modules to compute redistribution of the stress on the deformed RBC and subsequent re-deformation of the cell until an equilibrium state was reached [20], The Comsol multiphysics^TM software permits computing multiple physical problems in the same package and in a fully coupled manner. Our FEM analysis showed that the presence of the stress on the dimple part of the biconcave RBC surface could easily induce a large expansion of the RBC thickness. This unrealistic deformation was not observed in the experiments and can be prevented in the simulation only by reinforcing the constraints of constant RBC membrane volume and constant cytosol volume. Thus, we confirmed that in the optical tweezers experiments the volumes of the membrane and the cytosol keep constant. Comparison of the FEM predicted shape deformation with the experimental data provides the elastic parameters of the RBC.

2. Scattering of the RBC in dual-trap tweezers

We computed electromagnetic field of the twin trapping beams scattered by a biconcave RBC and the consequential radiation stress distribution on the RBC surface with the FEM. First, we modeled the vector field of the high numerical aperture (NA) trapping beams. The conventional Gaussian beam is a paraxial solution of the scalar Helmholtz equation and is not appropriate for describing a high NA laser beam. The vector Gaussian beam with high order corrections can describe accurately the tightly focused beam. The solution is described as a series of the vector spherical wave functions with the expansion coefficients computed with the point matching method in the T-matrix method [21]. However, a higher NA beam needs a higher order correction, but a higher order correction solution diverges at a shorter distance from the beam axis. For instance, a beam of NA = 1.25 requires the 7th order corrections, which converges within a distance of only 0.5λ from the beam axis [21] where λ = 1.06 μm is the wavelength of the beam. Therefore, this method would not be able to compute the field on the surface of a biological cell of 5-10λ size. In this paper we modeled the trapping beam as a spherical beam with the Gaussian intensity distribution in the cross sections [4]. The two trapping beams in the dual-trap optical tweezers propagating along the -z axis, shifted by a distance S along the x-axis and focused at (x = ± S/2, y = z = 0) have the electric field amplitude E₁ and E₂ expressed respectively as

E_{1, 2} (x, y, z) = {\begin{matrix} \sqrt{\frac{480 P}{n}} \frac{1}{W (z)} \exp [- \frac{{(x \pm S / 2)}^{2} + y^{2}}{W {(z)}^{2}}] \exp {j \frac{2 π}{λ} {[{(x \pm S / 2)}^{2} + y^{2} + z^{2}]}^{1 / 2}} (z > 0) \\ \sqrt{\frac{480 P}{n}} \frac{1}{W (z)} \exp [- \frac{{(x \pm S / 2)}^{2} + y^{2}}{W {(z)}^{2}}] \exp {- j \frac{2 π}{λ} {[{(x \pm S / 2)}^{2} + y^{2} + z^{2}]}^{1 / 2}} (z \leq 0) \end{matrix}

where

P = I_{0} (π W_{0}^{2}) / 2

is the beam power with the on-axis beam intensity I₀ and the beam waist radius

W_{0} = λ / π \tan (\sin^{- 1} (N A / n))

with the index of refraction n_1,2 of the media. The electrical field E is related to the beam intensity as

{| E |}^{2} = 2 I η

for the impedance

η = η_{0} / n

with η₀ = 120π Ω, the beam radius

W (z) = W_{0} {(1 + {(z / z_{0})}^{2})}^{1 / 2}

and the Rayleigh range

z_{0} = π W_{0}^{2} / λ

. The plus and minus signs of the phase of the spherical waves in the last term in Eq. (1) correspond to the fact that the beam is convergent when z > 0 and divergent when z ≤ 0. The polarization ofis in the plane normal to the wave vector of the spherical wavefront. Its components may be determined by the geometrical consideration. Note that in this vector field model the beam has a non-zero component along the beam axis (i.e. E_z ≠ 0). For instance, for a radial polarization beam the electrical field components for beams 1, and 2 can be expressed as

\begin{array}{l} E_{x 1, 2} = E_{1, 2} \cos ϕ_{1, 2} \sin θ_{1, 2} \\ E_{y 1, 2} = E_{1, 2} \cos ϕ_{1, 2} \cos θ_{1, 2} \\ E_{z 1, 2} = E_{1, 2} \sin ϕ_{1, 2} \end{array}

where ϕ is the zenith angle and θ the azimuth angle. The calculations reported in this paper were all done with the radial polarization beams. Other polarizations may be modeled in a similar manner.

The two trapping beams in the optical tweezers without the presence of the cell constitute the background field in the FEM. Then, a biconcave RBC is introduced into the background field to compute the scattering of the cell. Experimental observation showed that a stably trapped biconcave RBC was flipped such that its platelet surface was oriented parallel to the plane spanned by the axes of the two trapping beams [4], e. g. the x-z plane in Fig. 1(a), because in this orientation the cell, whose refractive index is higher than that of the surrounding medium, has the maximum volume within the laser beams, such that the potential energy associated to the light induced torque is minimum, where $\vec{p} = α \vec{E}$ is the polarization vector and α is the polarizability of the RBC, so that the trapping is stable.

Fig. 1 (a): Biconcave RBC in a dual-trap optical tweezers; (b) A quarter of the model used to compute light-scattering by the RBC.

Download Full Size | PDF

The diameter, thickness, surface area and volume of a healthy human RBC vary from sample to sample, following the normal distributions with a standard deviation of about 10%, and change with the osmotic pressure [22, 23]. In the simulation we used the formula for the biconcave shape given at the tonicity of 300 mOsm [2, 22]

y = \pm 0.5 R_{0} \sqrt{1 - \frac{(x^{2} + z^{2})}{R_{0}^{2}}} [c_{0} + c_{1} \frac{x^{2} + z^{2}}{R_{0}^{2}} + c_{2} \frac{{(x^{2} + z^{2})}^{2}}{R_{0}^{4}}]

where the radius of the RBC disk in the x-z plane R₀ is set to 3.90 µm and (c₀, c₁, c₂) = (0.207161, 2.002558, −1.122762).

In our model the geometrical center of the biconcave RBC was at the origin of the coordinate system and the focal points of the trapping beams are on the x-axis. Although a trapped cell can be pushed away from the focus of the trapping beam, we ignored this slight offset to simplify the numerical calculation. Ideally, in the numerical simulation one could leave the cell floating in the buffer and compute all the cell scattering, radiation stress, cell deformation and total force to determine the cell equilibrium position with the FEM. However, this will require a too large domain of calculation and a huge computer memory.

The FEM model to compute the cell scattering consisted of the cytosol and the membrane of the RBC, the buffer and the perfect matched layer (PML), as shown in Fig. 1(b). We set the same refractive index n₁ = 1.378 for RBC cytosol and membrane, and n₂ = 1.335 for the buffer solution at the wavelength λ = 1.06 μm and of the beam power 40 mW. The PML is an additional domain to absorb any radiation incident on the PML without producing reflections. We set a PML containing 5 spherical swept mesh layers with a total thickness of 0.6 μm, with the internal layer radius of 4.5 μm and the relative permittivity ε_r = 1.335². All the materials in the model have electric conductivity σ = 0 (S/m), and relative permeability μ_r = 1. The model was meshed with a total of 363,548 mesh elements in the volume of about 555 µm³ including the tetrahedral, triangular, prism, quadrilateral, edge and vertex elements. The mesh grid size was adaptively set in the range from λ/10 to λ/4. The total field as an addition of the background field and the scattered was computed by the Comsol RF module 4.2a^TM.

Figure 2(a) shows the background field on a biconcave RBC of diameter 7.80 μm in the dual-trap optical tweezers with the two beams separated by a distance S = 6.6 μm; the corresponding total electrical field is shown in Fig. 2(b). On the cell surface the background field magnitude has the highest value close to the extremities of the cell and the lowest value in the cell’s horizontal central sections. The background field in the cross section of x-y plane and x-z plane is shown in Fig. 3(a) and 3(c), respectively. Without taking into account the scattering by the cell the field is focused at x = ± S/2, y = 0, and z = 0, and is symmetric with respect to the x-y plane. When the cell scattering was taken into account, due to the dioptric power of the curved RBC/buffer interface, the trapping beams are partially scattered towards the cell’s center to form two weak spots at |x|< S/2 beside the initial focal spots at x = ± S/2, as perceptible in the cross sections of the cell in the x-y plane in Fig. 3(b) and in the x-z plane in Fig. 3(d). There is also a small shift of the maximum intensity spots within the focal volume of the beams in the + z-direction opposite to the beam propagation (shorter focal length), as can be seen in Fig. 3(d). As a result, the beams incident from the upper side of the cell are more dispersed after the focal points inside the cell when exiting from the lower side of the cell as shown in Fig. 2(b). In this simulation the absorption of the cytosol was not considered, as the absorption constant of the hemoglobin at the wavelength λ = 1.06 µm is very small (in the order of 10⁻⁶ [24]); we set the imaginary part of the complex refractive index of the RBC to be zero.

Fig. 2 (a) Background electrical field; (b) Total electrical field on an RBC of diameter 7.8 μm with two trapping beams separated by 6.6 μm.

Download Full Size | PDF

Fig. 3 Background field in the cross section of an RBC in the x-y plane (a) and in the x-z plane (c); Total field in the cross section of RBC in the x-y plane (b) and in the x-z plane (d). The RBC diameter = 7.8 μm, and the focal points of the two trapping beams were separated by 6.6 μm.

Download Full Size | PDF

3. Maxwell stress tensor

The radiation stress distribution at the cell/buffer interface was computed as the inner product of the Maxwell stress tensor and the normal of the interface and is expressed as

\vec{σ} = ({\overset{\leftrightarrow}{T}}_{2} - {\overset{\leftrightarrow}{T}}_{1}) \cdot \vec{n}

where the Maxwell stress tensor has the elements as [25]

T_{i j} = ε [E_{i} E_{j} + \frac{1}{ε μ_{0}} B_{i} B_{j} - \frac{1}{2} (E^{2} + \frac{1}{ε μ_{0}} B^{2}) δ_{i j}]

Considering the terms related to E only and

{(\overset{\leftrightarrow}{T} \cdot \vec{n})}_{i} = ε \sum_{j} E_{i} E_{j} n_{j} - (1 / 2) ε E^{2} n_{i}

we have

\vec{σ} = ({\overset{\leftrightarrow}{T}}_{2} - {\overset{\leftrightarrow}{T}}_{1}) \cdot \vec{n} = (ε_{2} E_{2 n} {\vec{E}}_{2} - ε_{1} E_{1 n} {\vec{E}}_{1}) - (1 / 2) (ε_{2} E_{2}^{2} - ε_{1} E_{1}^{2}) \vec{n}

with the first term in the right-hand side of Eq. (6) becomes

\begin{matrix} (ε_{2} E_{2 n} {\vec{E}}_{2} - ε_{1} E_{1 n} {\vec{E}}_{1}) = ε_{2} E_{2 n}^{2} \vec{n} - ε_{1} E_{1 n}^{2} \vec{n} + ε_{2} E_{2 n} {\vec{E}}_{2 t} - ε_{1} E_{1 n} {\vec{E}}_{1 t} \\ = (ε_{2} E_{2 n}^{2} - ε_{1} E_{1 n}^{2}) \vec{n} \end{matrix}

as according to boundary conditions

{\vec{E}}_{2 t} = {\vec{E}}_{1 t}

and

ε_{2} E_{2 n} = ε_{1} E_{1 n}

, so that

\begin{array}{l} \vec{σ} = [(ε_{2} E_{2 n}^{2} - ε_{1} E_{1 n}^{2}) - (1 / 2) (ε_{2} E_{2}^{2} - ε_{1} E_{1}^{2})] \vec{n} \\ = (1 / 2) [(ε_{2} E_{2 n}^{2} - ε_{1} E_{1 n}^{2}) - (ε_{2} E_{2 t}^{2} - ε_{1} E_{1 t}^{2})] \vec{n} \\ = (1 / 2) (ε_{1} - ε_{2}) [(\frac{ε_{1}}{ε_{2}}) E_{1 n}^{2} + E_{1 t}^{2}] \vec{n} \end{array}

where

E^{2} = E_{n}^{2} + E_{t}^{2}

and the boundary condition

E_{2 n} = E_{1 n} (ε_{1} / ε_{2})

were used. The two terms in Eq. (5) with the magnetic intensity are treated in the similar ways, resulting in zero, as

μ_{1} = μ_{2}

Eq. (7) is the same as that given in [26], except a factor of ½. According to Eq. (7) the radiation stress in the isotropic media is always normal to the interface and the stress can be computed from the normal and tangent components E_n and E_t of the electrical field in one medium only.

In the dual-trap optical tweezers the two trapping beams can interfere in their overlapping region, when the two beams coming from the same laser source. There is no such an interference when two individual lasers [5] or a single jumping laser beam [4] was used in the dual-trap optical tweezers. To avoid the beam interference in the simulation we computed the scattered fields of each single beam and the radiation stress distribution associated with each beam in two separated FEM modules, and then summed up two stress distributions on the cell surface.

Figure 4 shows the 3D stress distribution on the surface of a biconcave RBC in the dual-trap optical tweezers as a function of the beam separation distance S, varying from 3.1 to 7.3 μm. According to this result, there are two regimes of the stress distribution: for small beam separations the stress is mostly concentrated in the central dimple part of the RBC surface with high values, as shown in the first row in Fig. 4 for S ≤ 5.2 μm. This is conceivable because at side-incidence to the biconcave RBC the high NA trapping beams can actually reach the central pit of the RBC. In this regime, the dimple part is stretched outwards from the x-z plane, and there is slight shrinking of the RBC in the ± x directions, as observed in the experiments [4]. For large beam separations, the stress is mostly concentrated in the outer part of the biconcave surface. When increasing the beam separation the stress moves towards the two extreme ends of the cell in the ± x-axis and the maximum stress value increases. The larger the beam separation, the larger the maximum stress value and the closer the maximum stress spot to the x-axis, as shown in the second row in Fig. 4 for S ≥ 5.9 μm. Most experiments for stretching the RBC should work in this regime. Of cause, when the beam separation is larger than the diameter of the RBC disk the stress drops rapidly. The critical beam separation between the two stress regimes may be estimated approximately as the diameter of the circle of maximum thickness in the ± y-direction of the biconcave shape, given by 2r = 5.6 μm with $r = \sqrt{x^{2} + z^{2}}$ by setting the derivative dy/dr = 0 in Eq. (3).

Fig. 4 Stress distributions on a biconcave RBC trapped in the dual-trap optical tweezers. Beam separation S (in μm) is given under each figure. Value of the stress is referred to the color bar. Value of the stress was halved for S = 7.3 μm.

Download Full Size | PDF

Our 3D numerical result confirms the approximate estimation proposed in [4], which was based on the 2D geometrical optics with the biconcave RBC modeled as a disc platelet. However, the stress in the dimple part of the RBC was ignored in the 2D model [4].

4. Deformation of the RBC

Deformation of the RBC under the radiation stress distribution on the surface was computed by FEM with the continuum mechanics theory, which analyzes the kinematics and the interrelations among surface force, internal contact force, finite strain, displacement field and deformation in the continuum materials based on fundamental physical laws and constitutive relations. The RBC wall consists of a phospholipid bilayer and an underlying spectrin molecular network. The total wall can be mechanically modeled as an effective continuum membrane [2], which bears the external load and contains the liquid cytosol inside. The deformability of the membrane is an important characteristic of the RBC. We modeled the RBC membrane as a hyperelastic incompressible Mooney-Revlin material [27], defined by the following strain energy density function:

W = C_{10} ({\bar{I}}_{1} - 3) + C_{01} ({\bar{I}}_{2} - 3) + \frac{1}{2} κ {(J - 1)}^{2}

where, with the definition of the deformation gradient tensor F and of the left and right Cauchy-Green deformation tensor

C = F^{T} F

and

B = F F^{T}

,

{\bar{I}}_{1} = J^{- 2 / 3} I_{1}

and

{\bar{I}}_{2} = J^{- 3 / 4} I_{2}

is the first and second invariant of C respectively with

I_{1} = t r a c e (B) = λ_{1}^{2} + λ_{2}^{2} + λ_{3}^{2}

I_{2} = \frac{1}{2} [I_{1}^{2} - t r a c e (B^{2})] = λ_{1}^{2} λ_{2}^{2} + λ_{2}^{2} λ_{3}^{2} + λ_{3}^{2} λ_{1}^{2}

where λ_i (i = 1, 2, 3) is the principal extension ratio along the three principal axes and

J^{2} = \det (B) = λ_{1}^{2} λ_{2}^{2} λ_{3}^{2}

For incompressible material, the volume ratio J = 1. On the right-hand side of Eq. (8) the first term is related to the linear elastic deformation with the material constant C₁₀ set to 10 Pa, by comparing with the experimental measurement of the cell deformation. The second term is related to the non-linear elastic deformation. We set C₀₁ = (1/30) C₁₀ [28]. The third term, related to the bulk modulus κ, denotes the volumetric strain energy density associated with the volume deformation. To keep the membrane volume constant we set an artificially high value of κ = 10⁹ Pa, so that J approached 1 and the volume was constant. We did not set the constraint for the constant membrane surface area [29].

Our experiment on the deformation simulation showed importance of the constant membrane volume constraint. In the small trapping beam separation regime, where the radiation stress is more concentrated in the dimple part of the biconcave RBC surface with high values, as shown in Fig. 4 for S = 3.1, 3.8 and 4.5 μm, without the constant membrane volume constraint this stress can induce very large expansions of the cell thickness, which were not observed in the experiments [4]. This unrealistic deformation can be avoided in the simulation only when the Mooney-Rivlin solid model was used with high value of the bulk modulus κ. The equivalent continuum membrane is usually set to a thickness of 40 nm, e.g. 1% of the radius of the RBC disk [20]. Even though the phospholipid bilayer of the real RBC membrane is about 4-8 nm thick [23], to reduce numerical error in the FEM calculation of the membrane volume we increased the thickness of the membrane to h₀ = 200 nm to ensure even mesh grid and good mesh qualities within the membrane thickness and the cytosol. In the linear elastic deformation regime, the scaling of the membrane thickness did not affect the computed value of the membrane shear modulus μ₀ = G₀h₀. Increasing the value of h₀ would lead to less deformation. To produce the same amount of deformation as that measured in the experiment, the initial bulk shear modulus G₀ = 2C₁₀ should decrease in the same proportion, so that the membrane shear modulus μ₀ does not change.

The total volume enclosed by the cell membrane should be kept constant, as the osmotic exchange of material across the RBC membrane is unlikely during the short time of optical tweezers experiment. To keep the cell internal volume constant in the FEM simulation we filled the cell by the cytosol, which was set as the viscoelastic material with an extremely low Young’s modulus E = 1 Pa and Poisson ratio v = 0.49. We did not set the gravity force and did not set the time dependence in the properties of the cytosol and membrane for simulating the RBC deformation in an equilibrium state.

Although the biconcave shape and the two trapping beams are symmetrical with respect to the x-z plane and the y-z plane, the radiation stress distribution is not symmetrical to the x-y plane due to the cell scattering of the incident beams. As a result, the total force applied to the RBC was not zero and the whole cell can move in the space. To obtain the cell deformation in the FEM simulation we had to eliminate its movement by setting artificial constraints that the edges of the RBC in the three x-z, y-z and x-y planes can only move in the same planes, but cannot float freely in the space. These constraints can leave some errors in the predicted shape deformation of the cell.

5. Results

In the simulation the scattered electromagnetic field was computed in the Comsol RF modules^TM with the frequency domain solver at the laser beam frequency. The cell deformation was computed in the Comsol Structural Mechanics module^TM by the stationary solver. The two physics problems are coupled: once the cell’s morphology is deformed, the radiation stress distribution on the cell surface may change instantaneously. However, the reaction of the cell to the new stress distribution may take some time due to its viscosity. The subsequent cell deformation in turn will induce new stress redistribution. The process can be simulated step-by-step with Comsol Multiphysics^TM. The final equilibrium state of the cell may be obtained with the recursive loop as shown in Fig. 5. Started from the initial geometrical shape, we computed the scattered field, the stress distribution and cell deformation. Then, we updated the cell shape with the computed displacement vectors (u, v, w) to compute the next iteration. The iterations continued until convergence, when the residual error in each module was smaller than a prescribed tolerance value.

Fig. 5 Flow chart to compute RBC deformation in the equilibrium state with coupled RF and Continuum Mechanics Comsol^TM modules.

Download Full Size | PDF

Figure 6 shows the 3D deformed shape of an originally biconcave-shaped RBC in the equilibrium state computed by the FEM. In the stretching regime with the two beam separation S = 7.3 μm the RBC was elongated along the ± x axis from an initial value of 7.80 μm to a final value of 8.54 μm which is comparable with the experimental data [4]. As the hyperelastic coefficient C₁₀ = G₀/2 [2, 28], the setting value of C₁₀ = 10 Pa corresponds to the shear modulus of the membrane μ₀ = G₀h₀ = 4 μN/m. In the regime with small beam separation S = 3.80 μm, the radiation stress was concentrated in the dimple part of the biconcave shape, the RBC was not elongated but was slightly shrunken in the x-direction to 7.61 μm without large expansion in the ± y-direction. This is in agreement with the experimental observation [4].

Fig. 6 3D deformation of an RBC in dual-trap tweezers. First row: beam separation: S = 7.3 µm; Second row: S = 3.8 µm. Color bar: displacement positive (outward) or negative (inward).

Download Full Size | PDF

6. Conclusion

We have built a physical model with Comsol^TM FEM modules to simulate the deformation of a biconcave RBC trapped directly by two optical trapping beams. The light scattering by the RBC at the side-incidence was calculated with the solution of the Maxwell’s equations. The radiation stress distribution was computed with the Maxwell stress tensor from the total field. The morphological deformation of the cell at the equilibrium state was computed with the continuum mechanics module, which is fully coupled with the electromagnetic module in a recursive manner. The hyperelastic coefficients and the shear modulus of the RBC membrane were obtained through the comparison with the experimental deformation data. It is important to constraint the constant membrane volume and the constant cytosol volume in the simulation in order to obtain the simulated deformation results comparable with experiments. Our model can be also used to compute the viscoelasticity and other mechanical properties of the biological cells manipulated by the optical tweezers.

Acknowledgments

Yu and Sheng are supported by the Discovery research program of Natural Sciences and Engineering Research Council (NSERC) of Canada. Chiou is supported in part by the UST-UCSD International Center of Excellence in Advanced Bioengineering sponsored by the Taiwan National Science Council I-RiCE Program and by the Ministry of Education (The Top University Project).

References and links

1. 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]

2. A. Ashkin, J. M. Dziedzic, and T. Yamane, “Optical trapping and manipulation of single cells using infrared laser beams,” Nature 330(6150), 769–771 (1987). [CrossRef] [PubMed]

3. M. Khan, H. Soni, and A. K. Sood, “Optical tweezers for probing erythrocyte membrane deformability,” Appl. Phys. Lett. 95(23), 233703 (2009). [CrossRef]

4. 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]

5. M. Kinnunen, A. Kauppila, A. Karmenyan, and R. Myllylä, “Effect of the size and shape of a red blood cell on elastic light scattering properties at the single-cell level,” Biomed. Opt. Express 2(7), 1803–1814 (2011). [CrossRef] [PubMed]

6. A. C. De Luca, G. Rusciano, R. Ciancia, V. Martinelli, G. Pesce, B. Rotoli, L. Selvaggi, and A. Sasso, “Spectroscopical and mechanical characterization of normal and thalassemic red blood cells by Raman Tweezers,” Opt. Express 16(11), 7943–7957 (2008). [CrossRef] [PubMed]

7. E. V. Lyubin, M. D. Khokhlova, M. N. Skryabina, and A. A. Fedyanin, “Cellular viscoelasticity probed by active rheology in optical tweezers,” J. Biomed. Opt. 17(10), 101510 (2012). [CrossRef] [PubMed]

8. Y. Z. Yoon, J. Kotar, A. T. Brown, and P. Cicuta, “Red blood cell dynamics: from spontaneous fluctuations to non-linear response,” Soft Matter 7(5), 2042–2051 (2011). [CrossRef]

9. I. Sraj, A. C. Szatmary, D. W. M. Marr, and C. D. Eggleton, “Dynamic ray tracing for modeling optical cell manipulation,” Opt. Express 18(16), 16702–16714 (2010). [CrossRef] [PubMed]

10. J. H. Zhou, M. C. Zhong, Z. Q. Wang, and Y. M. Li, “Calculation of optical forces on an ellipsoid using vectorial ray tracing method,” Opt. Express 20(14), 14928–14937 (2012). [CrossRef] [PubMed]

11. J. P. Mills, M. Diez-Silva, D. J. Quinn, M. Dao, M. J. Lang, K. S. W. Tan, C. T. Lim, G. Milon, P. H. David, O. Mercereau-Puijalon, S. Bonnefoy, and S. Suresh, “Effect of plasmodial RESA protein on deformability of human red blood cells harboring Plasmodium falciparum,” Proc. Natl. Acad. Sci. U.S.A. 104(22), 9213–9217 (2007). [CrossRef] [PubMed]

12. G. J. C. G. M. Bosman, “Erythrocyte aging in sickle cell disease,” Cell. Mol. Biol. (Noisy-le-grand) 50(1), 81–86 (2004). [PubMed]

13. Y. K. Park, M. Diez-Silva, G. Popescu, G. Lykotrafitis, W. S. Choi, M. S. Feld, and S. Suresh, “Refractive index maps and membrane dynamics of human red blood cells parasitized by Plasmodium falciparum,” Proc. Natl. Acad. Sci. U.S.A. 105(37), 13730–13735 (2008). [CrossRef] [PubMed]

14. A. Krantz, “Red cell-mediated therapy: Opportunities and challenges,” Blood Cells Mol. Dis. 23(1), 58–68 (1997). [CrossRef] [PubMed]

15. A. Karlsson, J. He, J. Swartling, and S. Andersson-Engels, “Numerical simulations of light scattering by red blood cells,” IEEE Trans. Biomed. Eng. 52(1), 13–18 (2005). [CrossRef] [PubMed]

16. J. Lim, H. Ding, M. Mir, R. Zhu, K. Tangella, and G. Popescu, “Born approximation model for light scattering by red blood cells,” Biomed. Opt. Express 2(10), 2784–2791 (2011). [CrossRef] [PubMed]

17. T. Wriedta, J. Hellmers, E. Ereminab, and R. Schuh, “Light scattering by single erythrocyte: Comparison of different Methods,” J. Quant. Spectrosc. Ra. 100(1-3), 444–456 (2006). [CrossRef]

18. J. Q. Lu, P. Yang, and X. H. Hu, “Simulations of light scattering from a biconcave red blood cell using the finite-difference time-domain method,” J. Biomed. Opt. 10(2), 024022 (2005). [CrossRef] [PubMed]

19. J. T. Yu, J. Y. Chen, Z. F. Lin, L. Xu, P. N. Wang, and M. Gu, “Surface stress on the erythrocyte under laser irradiation with finite-difference time-domain calculation,” J. Biomed. Opt. 10(6), 064013 (2005). [CrossRef] [PubMed]

20. S. Rancourt-Grenier, M. T. Wei, J. J. Bai, A. Chiou, P. P. Bareil, P. L. Duval, and Y. Sheng, “Dynamic deformation of red blood cell in dual-trap optical tweezers,” Opt. Express 18(10), 10462–10472 (2010). [CrossRef] [PubMed]

21. P. B. Bareil and Y. Sheng, “Modeling highly focused laser beam in optical tweezers with the vector Gaussian beam in the T-matrix method,” J. Opt. Soc. Am. A 30(1), 1–6 (2013). [CrossRef] [PubMed]

22. E. Evans and Y. C. Fung, “Improved measurements of the erythrocyte geometry,” Microvasc. Res. 4(4), 335–347 (1972). [CrossRef] [PubMed]

23. X. C. Fung, Biomechanics: Mechanical Properties of Living Tissue, 2nd ed. (Springer, 1993).

24. M. Friebel and M. Meinke, “Determination of the complex refractive index of highly concentrated hemoglobin solutions using transmittance and reflectance measurements,” J. Biomed. Opt. 10(6), 064019 (2005). [CrossRef] [PubMed]

25. J. D. Jackson, Classical Electrodynamics, 3rd ed. (Wiley, 1998).

26. K. Okamoto and S. Kawata, “Radiation force exerted on subwavelength particles near a nanoaperture,” Phys. Rev. Lett. 83(22), 4534–4537 (1999). [CrossRef]

27. E. A. Evans, “A new material concept for the red cell membrane,” Biophys. J. 13(9), 926–940 (1973). [CrossRef] [PubMed]

28. C. Y. Chee, H. P. Lee, and C. Lu, “Using 3D fluid-structure interaction model to analyze the biomechanical properties of erythrocyte,” Phys. Lett. A 372(9), 1357–1362 (2008). [CrossRef]

29. R. Skalak, A. Tozeren, R. P. Zarda, and S. Chien, “Strain energy function of red blood cell membranes,” Biophys. J. 13(3), 245–264 (1973). [CrossRef] [PubMed]

Three-dimensional light-scattering and deformation of individual biconcave human blood cells in optical tweezers

Abstract

1. Introduction

2. Scattering of the RBC in dual-trap tweezers

3. Maxwell stress tensor

4. Deformation of the RBC

5. Results

6. Conclusion

Acknowledgments

References and links

Cited By

Figures (6)

Equations (13)

Optics Express