Abstract
Magnetic twisting cytometry probes mechanical properties of an adherent cell by applying a torque to a magnetic bead that is tightly bound to the cell surface. Here we have used a threedimensional finite element model of cell deformation to compute the relationships between the applied torque and resulting bead rotation and lateral bead translation. From the analysis, we computed two coefficients that allow the cell elastic modulus to be estimated from measurements of either bead rotation or lateral bead translation, respectively, if the degree of bead embedding and the cell height are known. Although computed strains in proximity of the bead can be large, the relationships between applied torque and bead rotation or translation remain virtually linear up to bead rotations of 15°, above which geometrical nonlinearities become significant. This appreciable linear range stands in contrast to the intrinsically nonlinear forcedisplacement relationship that is observed when cells are indented during atomic force microscopy. Finally, these computations support the idea that adhesive forces are sufficient to keep the bead firmly attached to the cell surface throughout the range of working torques.
 elastic modulus
 cell stiffness
 bead embedding
 cell height
 magnetic twisting cytometry
cell deformability plays an important role in a variety of cell functions, including contraction, spreading, crawling, and wound healing (46, 13). To measure cell deformability and its changes, our laboratory has made extensive use of magnetic twisting cytometry (5, 6, 14, 15, 20, 22, 23, 29, 35, 3739). Magnetic twisting cytometry probes mechanical properties of an adherent cell by applying a torque to a magnetic bead that is tightly bound to the cell surface. Most studies that use this method have reported the relationship between the applied torque and resulting angular rotations. More recent studies have reported the relationship between the applied torque and resulting lateral bead displacements (10,11). With either approach, data analysis reported thus far has emphasized relative changes in cell stiffness, or an “apparent stiffness,” rather than the absolute elastic modulus.
In this report, we have used threedimensional finite element models to compute the relationships between the applied torque and resulting cell deformation, bead rotation, and lateral bead translation. We have assessed the effects of different degrees of bead embedding and cell height within a geometrically linear range of cell deformation. From these relationships, the elastic modulus of the cell can be computed if the beadcell geometry is known. Assuming a linear elastic material and using geometrically nonlinear finite element analysis, we also assessed the limit up to which the linear analysis can be used.
METHODS
Finite Element Model, Boundary Conditions, and Model Parameters
Geometry.
The cell was modeled as a threedimensional slab with locally constant height and a lateral extent of 50 bead diameters (Fig.1 B). Justification for the assumption of a locally constant cell height is addressed in thediscussion. We also considered a spherical bead embedded in an infinite elastic space that could be compared with an exact analytic solution (28).
Finite element formulation.
The finite element program PAK (17, 19) was used to calculate the rotation and displacement of a rigid sphere partially or completely embedded in the deformable continuum (cell). The deformable continuum is defined as elastic solid of finite or semiinfinite thickness, with the bottom surface fixed to the rigid surface. For a torque applied along the xaxis (Fig.1 A), we calculated displacements, strains, and stresses within the material (i.e., cell) and also the bead translation and rotation. We assumed no slipping between the bead surface and the cell body at the interface.
The essence of the finite element method consists of division of a continuum into subdomains, called finite elements, in which an approximate displacement field u is defined as
We assumed that bead stiffness is much larger than cell stiffness, in which case, considering bead as a rigid body, the bead motion is defined by displacement of its center and its rotation. We impose the condition that displacements at the interface between the rigid bead (u
_{bead}) and the deformable cell (u
_{cell}) are the same
In the finite element analysis, we used threedimensional, solid, isoparametric, linear eightnode elements. Typically, the total number of nodes was ∼4,500 with ∼15,000 elements.
We performed two analyses; in both we assumed that the material was linearly elastic. The first analysis was geometrically linear in the sense that it did not take into account geometrical distortions induced by bead motion. The second analysis was geometrically nonlinear in that it did take into account such distortions. To do this latter analysis, we included quadratic terms in the strain calculations and change of geometry. Thus the Cauchy (true) stresses were related to the Almansi strains in the linear constitutive law.
We used a Poisson ratio, υ = 0.49, a value that does not cause pressure locking of the elements. We also performed calculations with υ = 0.45 and 0.499 and found only minor changes in sphere rotation and displacement (differences of <0.5%). The stress field corresponding to υ = 0.49 was altered by only a few percent in the case when υ was equal to 0.45. Thus our calculations with υ = 0.49 are reasonably accurate.
The finite element mesh.
We point out that the stress everywhere on the beadcell interface is expected to be bounded. This result is nontrivial, biologically important, and impacts the question of what is a sufficiently fine mesh for reasonable finite element approximations. The stress is bounded, even at the extremes of the beadcell interface, because the sphere acts on the deformable cell as a rigid body without sharp edges (24). This is in contrast to the case of a rigid punch with sharp edges, which causes infinite stresses in elastic deformable solid in the proximity of the punch edges.
In the case of a bead embedded by only 10% of its diameter, we performed calculations with two meshes, one of which had 50% more elements than the other. In regions of maximum stress, we found differences approaching ∼2.5%, with much smaller discrepancies elsewhere. The corresponding differences in computed bead displacements and rotations were <0.1%. Thus the use of the finer mesh was not necessary.
The dependency of the results on the mesh size in the regions of high deformations is important, however, in analysis of geometrical nonlinearities. We observed negative Jacobians (i.e., the determinant of the transformation matrix between Cartesian and isoparametric coordinates) in a few of the most distorted elements when bead rotation exceeded 15°. An extreme case, in which Jacobians are still positive but small, is shown in Fig. 1 B. The contribution of these distorted elements to the overall response was negligible, however. In the linear analysis, the Jacobians remained constant and are independent of the element distortion; in that case, we limited our calculations to specific torques (T_{s}) that were 25% smaller than those used in the nonlinear analysis. Numerical computations for the case of a bead embedded in an infinite medium (Fig. 2 A) agreed within 0.5% with the corresponding analytic solution of PhanThien (28).
Parameter values.
The diameter of the ferromagnetic bead was taken as 4.5 μm. The elastic material parameters of the deformable continuum (i.e., the cell) were taken as shear modulus (G) = 1,000 Pa and υ = 0.49; the G is related to the Young's modulus (E) by the relationship G = E/[2(1 + υ)]. We studied beads embedded to a depth of 5, 10, 25, 75, and 100% of the bead diameter, and cell heights of 1 and 5 μm and semiinfinite.
Two different definitions for T_{s} have been used in the literature, and they differ by geometric scaling factor (6 for spherical beads) (10, 11, 20, 36, 38). Here we follow Fabry et al. (10, 11) by defining T_{s} as applied magnetic torque (T_{m}) per unit of bead volume. For the linear finite element analysis, the T_{s} acting on the bead is taken to be T_{s} = 60 Pa.
Calculation of Apparent Cell Stiffness from Measured Rotation or Lateral Translation
The magnetic twisting field (i.e., the applied magnetic moment) causes both a rotation and a lateral translation of the beads (Fig. 1). The bead rotation is φ, and the normalized lateral translation isd* = d/R [i.e., the bead lateral translation, d, expressed as a fraction of the bead radius (R)]. We define apparent cell stiffness as G_{φ} = T_{s}/(κφ) obtained from the bead rotation, or apparent cell stiffness as G_{d} = T_{s}/(κd*) from the bead lateral translation. Here κ is a shape factor (6 for spherical beads). We will show below that G_{φ} and G_{d} strongly depend on the degree of bead internalization and to a lesser degree on the cell height.
We calculated the relationship between the cell G and the apparent stiffnesses G_{φ} and G_{d} from torquerotation and torquedisplacement dependencies obtained in our analysis. We defined coefficients α and β by the relationships G_{φ} = αG and G_{d} = βG and then computed α and β each as a function of the depth of bead embedding and cell height.
RESULTS
Deformation, Stress, and Traction
Cell deformation and the stress fields in the plane perpendicular to the magnetic torque passing through the center of the bead are shown for four representative cases (Fig. 2). In all cases, the bead is taken out of the view, but the stresses at the cellbead interface are shown. We show the equivalent (von Mises) stress (1), which is the overall measure of maximum shear stress. The simplest case is rotation about the xaxis of a bead in an infinite medium (Fig. 2 A). In this case, there is no bead translation: the strain and stress fields and the interfacial tractions are uniform along the circumference in the plane perpendicular to the T_{m}. Stress decays rapidly in the radial direction (i.e., ∼R/r ^{3}), from the maximum value at cellmicrobead interface to ∼15% of that value at the radial distance of ∼2 μm from the microbead surface. Stress at the beadcell interface (shown inside the white circle in Fig.2 A) decays to zero at the central point that lies on the spherical surface and the center of the axis of rotation. The maximum value of the interfacial traction is in the plane passing through the center of the bead; for a bead twisted by a T_{s} of 60 Pa, the traction is equal to 30 Pa.
In contrast to the case of a bead in an infinite space, for a bead 50% embedded in a cell 5 μm in height (Fig. 2 B), the stress distribution is highly nonuniform around the bead interface and is maximal near the intersection of the bead and the free apical surface of the cell. This stress nonuniformity is far more prominent when the bead is only 10% embedded (Fig. 2, C and D, for two different cell heights). For 50% embedding, the peak stress is about two times larger than that of the bead immerged in an infinite medium, whereas for 10% embedding, the peak stress is >25 times larger. Moreover, reduction of cell height from 5 to 1 μm (Fig. 2,C and D) further increases the peak stresses by ∼20%. The maximum stress is at the bead surface and coincides with maximum traction at the cellbead interface. The maximum tractions also increase rapidly for less embedded beads. For example, for 10% embedding, the maximum traction is ∼2.5 kPa (Fig.3). The effect of these local peak stresses and peak tractions at the beadcell interface is discussed below.
Cell Shear Modulus
Figure 4 A shows the relationship between the coefficient α and the degree of bead embedding for three different cell heights. With less bead embedding, α is much smaller than 1. For example, for 10% embedding, G_{φ} is only 5% of G (Fig. 4 A), but α increases toward unity as bead embedding increases. A decrease of the cell height below onehalf of the microbead diameter at fixed bead embedding increases α more significantly, because the smaller bead rotation is caused by greatly reduced cell deformation in proximity of the rigid substrate.
Figure 4 B shows the relationship between the coefficient β and the degree of bead embedding for three different cell heights. Interestingly, the values of β are only slightly larger than α for the small degrees of bead embedding, but they become progressively larger than α as bead embedding increases. For a large degree of embedding (>70%), β becomes much larger than 1 and unbounded for the bead embedded in an infinite medium, whereas α never exceeds 1. These differences in α vs. β for the same degree of microbead embedding follows from the fact that, for the completely immerged microbead, the rotation remains finite, whereas the lateral displacement approaches zero.
The coefficients α and β decrease progressively when the cell height decreases below one bead diameter but change relatively little when the cell height is larger than one bead diameter (Fig. 4).
Nonlinearity
Partially embedded beads show stress concentrations and associated large local deformations (Fig. 2, C and D). For example, the maximum strain for 10% embedded bead may reach 30% at T_{s} = 60 Pa and up to 50% at T_{s} = 80 Pa (Fig. 5). These maximum strain values are much larger than the usual small deformation limit of 4% (1,18). Moreover, these large deformations create geometrical nonlinearities. To our surprise, the nonlinear analysis (performed for 10% embedded bead) showed that the relationship between the T_{s} and the bead rotation and/or bead displacements is insensitive to these local geometrical nonlinearities and remains virtually linear up to a rotation angle φ = 15°. For φ > 15°, geometrical nonlinearities became significant.
DISCUSSION
Magnetic twisting cytometry probes mechanical properties of an adherent cell by applying a torque to a magnetic microbead that is tightly bound to the cell surface. Here we computed coefficients α and β that allow the cell elastic modulus to be estimated from measurements of either bead rotation or lateral bead translation, respectively, if the applied torque, the degree of bead embedding, and the cell height are known. In the discussion, we deal first with three important limitations of the computational model and then go on to discuss the implication of our results.
First, this theoretical analysis of cell deformation is based on a highly simplified geometry of a spherical bead partially embedded in a homogeneous isotropic linear elastic slab whose height is locally uniform. Real cell geometry is, of course, far more complex (40), with cell height varying along and across the adherent cell. The analysis shows, however, that fields of stress and strain induced by the imposed torque decay rapidly in the radial direction, varying roughly as the inverse cube of distance from the bead; within one bead diameter, the stress and strain become quite small (Figs. 2 and 5). This highly localized response implies that geometrical details at lateral distances of more than one bead diameter would be expected to have at most only a slight influence on the computed results. Therefore, the analysis would be expected to be valid for regions in which changes of cell height are small over distances comparable to one bead diameter.
Second, to perform these calculations, we used only one value for the G, a realistic one, but the results are applicable, nonetheless, for any value of the G. In the linear case, the assumed value of the G acts merely as a simple scale factor, setting the level of the stress field relative to that of the strain field. For a fixed applied torque, the stresses and tractions depend on the particular numerical value of the modulus. Note that spatial distribution of strains depends neither on the applied torque nor on the modulus. Much the same is true in the case of induced nonlinear distortions: the angular rotations required to induce a geometric nonlinearity would be unchanged, but the level of the torque required to induce that rotation would change.
Third, the computational analysis reported here was limited to static conditions and does not explicitly take into account dynamic properties of the cytoskeletal lattice, which could cause both frequency dependence of the elastic modulus and an appreciable loss modulus corresponding to internal frictional stresses (10). Before, we defined G_{φ}= αG and G_{d} = βG, and now we define the complex moduli G̃_{φ}(f),G̃_{d}(f), andG̃(f), where f is the frequency of the excitation. These moduli denote complex variables possessing real (inphase or elastic) components and imaginary (outofphase or frictional) components. Assuming that the cytoskeletal material is linear and incompressible, the computed results can be generalized such that G̃_{φ}(f) = αG̃_{φ}(f) andG̃_{d}(f) = βG̃_{d}(f). The coefficients α and β are independent of frequency, and they are identical to the coefficients shown in Fig. 4. It can be shown that these coefficients equally apply to the cell constitutive law that obeys linear structural damping (810, 12) or linear viscoelasticity. The analysis should not be expected to extend to the case of plastic responses, however.
The key relationships for coefficients α and β are shown in Fig. 4. These relationships are counterintuitive in an important regard that deserves special comment. The results shown in Fig. 4 are surprisingly insensitive to cell height. For example, the results for infinite cell height vs. those for 5μm cell height show only negligible differences for a small extent of bead embedding and, at most, only modest differences for 100% bead embedding. In the latter case (100% embedding of a 4.5μm bead in a 5μm thick slab), there is a rigid substrate with a noslip boundary condition only 0.5 μm below the bead. It might have been thought that such close proximity of the rigid substrate would severely constrain bead motion and to an extent that is far greater than the analysis actually predicts. The analysis shows, however, that the dominant parts of the normal and tangential tractions that balance the torque applied to the bead are generated well away from the regions that most closely approach the rigid substrate (Fig.3), thus resolving the paradox. However, the effect of the proximity of the rigid substrate becomes important when cell height is below onehalf of bead diameter.
The coefficients α and β (Fig. 4) allow calculation of the cell elastic modulus from either torqueangle or torquedisplacement measurements. The elastic modulus obtained in that way is comparable to the modulus obtained by other techniques. For example, the apparent cell stiffness of 60–100 Pa observed in smooth muscle cells (20) divided by the coefficient α = 0.1 (that corresponds to average bead embedding of 20%) gives a cell G in the range of 600–1,000 Pa, which is comparable to measurements of the G of fibroblasts by micropipette manipulation (34) and cell poking (27). A similar value of the G is observed in myocytes obtained by atomic force microscopy (32). Similarly, an apparent shear cell stiffness of endothelial cells of ∼10 Pa (38) divided by the same coefficient (0.1) becomes close to the G obtained from the micropipette aspiration measurements of ∼100 Pa (30, 31).
The analysis shows that the degree of bead embedding and the cell height profoundly influence bead rotation and translation. Therefore, observed variations of bead embedding, cell height variation with position of the bead on each cell (i.e., close to nucleus, where the smooth muscle cell height is >5 μm, and at the cell periphery where the cell height of the wellspread cell is <2 μm), may be the major factors contributing to the variability of the apparent stiffnesses calculated from the observed displacements of the individual bead displacements (10, 11). However, this variability can also be caused, at least in part, by differences in bead attachments to the cytoskeleton, by local stiffness differences within a cell, or differences between cells (different levels of activation, for example) that collectively contribute to the observed heterogeneity of apparent stiffness.
Our analysis gives a basis for an estimation of the range of linearity of the torquerotation relationship over a wide range of bead embedding. In continuum mechanics of linear elastic solids, large deformations usually induce significant nonlinearity in the stressstrain relation for strains >4% (1, 18). Note that a strain of 4% corresponds to a bead rotation in an infinite medium of only 1.5°. Surprisingly, we found that, in the case of 10% embedding, the torquerotation relationship remains linear up to a bead rotation of 15° or lateral translation of 20% of the bead radius, which, for a radius of 2.25 μm, yields a lateral translation of ∼450 nm. This finding is consistent with experimentally observed linear rheological behavior among smooth muscle cells (10,11). Although the maximum local strains approach 50%, which greatly exceeds the ordinary linearity limit of 4%, the region of these large strains is quite localized to the edges of the cellbead interface (Fig. 5), and, as such, the overall system response remains virtually linear. With an increase of degree of bead embedding, these regions of large strains are less localized, and the linear torquerotation limit must be gradually decreased (from bead rotation of 15° at 10% embedding to bead rotation of 1.5° in an infinite medium), regardless of a large decrease in the peak strains and stresses (Fig. 2). Despite this decrease of the linear limit, the maximum bead rotation of 15° is a sufficient criterion to ensure that a population of beads embedded >10% is within the linear torquerotation range.
The existence of an appreciable linear range using this method stands in contrast to the forcedisplacement relationship measured with atomic force microscopy, which exhibits no linear range due to intrinsic geometric distortions induced by contact of the probe tip with the cell (7).
It is important to note that the bead lateral translations become extremely small for largebead embedding (Fig.6). For example, in the limit of the resolution of optical readings of ∼5 nm (10, 11), the displacement (of 4.5μm bead) can be accurately measured for beads embedded <75%. Moreover, determination of dynamic moduli from torquedisplacement loops may further be limited by the optical resolution to beads embedded <50%.
Finally, in the literature, the interpretation of data obtained by magnetic twisting cytometry has rested on the assumption that receptorligand bonds have sufficient adhesion strength and numbers to keep the bead attached to the cell over the whole interface and ensure a noslip boundary condition. Based on the tractions computed here, is that a reasonable assumption? Typical receptorligand bond strength is 10 pN, and a typical receptorligand bond density has a range from several hundred up to 10,000 bonds/μm^{2} (3, 16, 25,26, 41, 42). To support induced translation due to bead twisting, the minimum receptor density for the completely embedded bead (Fig. 2 A) is only 3 bonds/μm^{2}, whereas, for a bead emerged 10%, the minimum required bond density may reach >250 bonds/μm^{2}. Therefore, the ligand adhesion strength may be sufficient to keep the bead firmly attached to the cell and the cell to the substrate. Consequently, peeling of the bead from cell surface or the cell from substrate would not be expected to occur (10, 11,20).
Conclusions
We have shown that the elastic modulus can be obtained from the apparent cell stiffness if the degree of bead embedding and cell height is known. The elastic modulus estimated in this way is consistent with estimates of the elastic modulus assessed by other techniques in muscle and nonmuscle cells (20, 33, 43). We have calculated the limiting maximum angle of rotation in terms of the degree of the bead embedding, for which the torqueangle of rotation and torquedisplacement relationships are linear. Finally, based on estimated bond densities and maximum T_{s} values that are used in these techniques, we found that the adhesion strength of integrin bonds seems to be sufficient to keep the bead firmly attached to the cell.
Acknowledgments
We gratefully acknowledge Dr. Ning Wang and Dr. James Butler for helpful discussions and critical review of this article.
Footnotes

This work was supported by National Heart, Lung, and Blood Institute Grants HL33009 and HL59682.

Address for reprint requests and other correspondence: S. M. Mijailovich, Physiology Program, Dept. of Environmental Health, Bldg. I, Rm. 1306B, Harvard School of Public Health, 665 Huntington Ave., Boston, MA 02115 (Email:smijailo{at}hsph.harvard.edu).

The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked “advertisement” in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

June 21, 2002;10.1152/japplphysiol.00255.2002
 Copyright © 2002 the American Physiological Society