The optical diagnostic technique is potentially capable of rapid, noninvasive assessment of biological cells and tissues at a cellular level. The basis for all these applications is a fundamental understanding of the interaction of light with biological materials, i.e., light propagation (including absorption and scattering) in biological cells and tissues. Light propagation in biological cells is dependent on the tissue’s morphological and biochemical structure. This kind of dependence provides the possibility of obtaining information about the biological medium from the measured scattered light.
Cytometry is a technique that makes use of the measured scattered light to distinguish cells and has been applied in cell diagnostics for many years. Physically, light scattering occurs when the refractive index (or dielectric constant) is different between two different media. Therefore, the scattered light depends only on the shape of the interface and the electric properties of each medium. Although a biological cell is complex in its functions for a biologist, the cell is thought of by a physicist as a continuum containing different organelles and subcomponents of organelles, each with a different refractive index. Although cytometry is a well-established research and clinical diagnostics tool that has seen many refinements, some issues remain unclear. One example is the selection of the optimal angle for measuring side scattered light. Generally speaking, the cell size dominates forward scattering and the organelles dominate the side scattering. In previous analyses, cells are modeled as regular uniform spheres. Therefore, these treatments are unable to predict the behavior of side scattering, and users of cytometry select the angle of side scattering based on experience. To address the determination of an optimal observation angle in cytometry, a complete consideration of light scattering from cells in realistic situations is necessary.
In the analysis of the propagation and scattering of light in biological cells, the basic equations that govern light propagation and scattering are Maxwell’s equations. A direct analytical solution for Maxwell’s equation for light scattering from cells in realistic situations is impossible. In cell diagnostics, knowledge of two aspects are needed: for a given cell shape and refractive index distribution, researchers hope to exactly predict the scattered light; and from measured scattering data, one hopes to deduce the cell shape and its internal structure. There are no simple mathematical models that link the expected scattered light to realistic biological cells—only a few empirical expressions exist.1 Some approximate theories are used (see Ref. 2 for examples and references therein), but their applications are limited.
With advances in computational techniques, the numerical treatment of Maxwell’s equations provides an accurate way to address the issue of light scattering. There are several numerical studies reported.3 4 5 The finite-difference time-domain (FDTD) method has been used to simulate light scattering in biological tissues and cells over a given angular range.4 The FDTD method was first developed for calculating radar cross sections and can be applied to the case of complex scatters.6 Typical biological tissues and cells have irregular shapes and complex structures. Therefore, the FDTD method is a good choice for analyzing light scattering in biological cells.
The motivation for the work is the development of a computational tool that can be used for a complete and accurate understanding of light scattering in biological cells. A 3-D FDTD code is developed to simulate the light scattering from biological cells in realistic situations for all angles. The code that is described can serve as a tool to predict the scattered light from arbitrary cells. Furthermore, the features of the scattered light pattern are studied in detail for particular cell types. These can advance our knowledge of options that may be available in two-parameter cytometry. A method for determining the optimal observation angles for side scattered light is suggested. The remainder of this work is arranged as follows: in Sec. 2, a complete physics description for the code is given; and in Sec. 3, features of the scattered light are discussed in the view of cell diagnostics in two-parameter cytometry. A method for selecting the best angle for side scattered light is suggested. This is followed by a conclusion.
Basic Equations and Methods
The FDTD method dates back to the 1960’s, when Yee first suggested his particular discrete scheme for Maxwell’s equations.7 In the last two decades, the FDTD method has rapidly progressed with developments in computational technology.6 In this section, we outline our implementation of the code. Figure 1 gives the geometry used in our FDTD simulations. The entire calculation domain consists of the boundary, the near-to-far (NTF) field interface and the total field/scattered field (TF/SF) interface. The scatters (i.e., tissues or cells) are assumed to be located inside the TF/SF interface. The cell structure is modeled by assigning values of the permittivity to the various cell components.
Maxwell’s Equations and Their Discretization
The propagation, absorption, and scattering of laser light in any media is completely described by Maxwell’s equations. In the absence of the absorption and additional sources, the equations read:7 the six components of the electric and magnetic fields are separated from one another by a distance equal to 1/2 times the grid spacing, as displayed in Fig. 2. The Yee algorithm is a second-order, explicit, center-difference scheme. If we select the computational domain such that only components of the electric field lie on the domain surfaces and the grid is numbered from 1 to ns (s=x,y,z), then the discrete form of the curl equations for the electric and magnetic fields is as follows: εi,j,k and μi,j,k are values of permittivity and permeability at point (i,j,k), respectively. Δt is the time step. Subscripts i, j, and k represent the grid numbers in x, y, and z directions. Δx, Δy, and Δz are the grid spacings in x, y, and z directions, respectively. Other components of the electric and magnetic fields may be obtained in a similar fashion. It has been proven that in the Yee discretization scheme, the divergence equations of Maxwell’s equations are automatically satisfied.6 The Courant condition for Eq. (3) is given by:6 6 that the spatial grid spacing must be less than λ/10, where λ is the wavelength of the incident light.
Although light actually propagates in an unbounded space, our computations can only be conducted in a limited spatial domain. Special boundary conditions must be applied to minimize light reflection from the boundary in order that the light appears to be propagating as if it were in an unbounded space. There have been numerous boundary conditions reported to address this issue. Three of the most frequently used boundary conditions are the Mur boundary condition,8 the Liao boundary condition,9 and the perfect matched layer (PML) boundary condition.10 In our code development we have tested all three of these boundary conditions. In the simulation presented in this work, the Liao boundary condition has been used.
Total Field/Scattering Field Formulation
The total field/scattering field (TF/SF) technique is used to input the incident light and to extract the scattering field from the total field. This technique simplifies the implementation of the scattered field calculations. The total field is calculated inside the TF/SF interface in Fig. 1. Outside the TF/SF interface, only the scattered field is calculated. From the linearity of Maxwell’s equations, we have:
Incident Light Calculation
One can calculate the incident light on the TF/SF interface either directly from the waveform formula, or by following the approach of a table look-up procedure introduced by Taflove,6 and Oguz and Gurel.11 We follow the second approach, which is better in reducing errors caused by numerical dispersion. When the incident light is a plane wave propagating along the positive z direction, values of the incident electromagnetic field are obtained from the 1-D Maxwell equations. The x or y components of the incident electromagnetic fields at each grid of the TF/SF interface are then easily obtained by using the polarization property of the incident light.
To avoid high-frequency components that may be introduced by a sudden onset of the incident wave, the Hamming function is used for describing the magnitude of the startup of the incident light:ΔT denotes the startup time. In our calculation, ΔT is selected to be ω/2, where ω is the angular frequency of the incident light.
Near-to-Far Field Transformation
In most physical situations, the interest is in the scattered light far from the scatterer itself. However, FDTD computations can only be used to obtain the scattered electric field and magnetic field close to the scatterer. The technique of equivalent electric/magnetic currents is used to calculate the scattered field at a distance by using the near-field data. The NTF interface is an enclosure surface that is virtually located in the scattered field domain, i.e., between the calculated domain boundary and the TF/SF interface. The radiation source inside the NTF interface can be equivalently represented by virtual electric currents J&vec;s and the magnetic currents M&vec;s, which flow on the NTF interface. These are defined by:E&vec; and H&vec; are the calculated electromagnetic fields on the NTF interface, and s&vec;0 is the unit vector normal to the interface. There are two strategies that can be used to implement the NTF field transformation: the time-domain approach and frequency-domain approach. We employ the frequency-domain NTF transformation. In this method, all components of electromagnetic fields are transformed from a time-domain value into a frequency-domain value by using the discrete Fourier transform (DFT). By using a spherical coordinate system (r,θ,ϕ), where θ is the angle between r&vec; and z&vec;, the far-field calculations are performed as follows: N&vec; and L&vec; are defined as: r′ and the field point r is given by ϕ, which can be approximated as: r′ cos ϕ=x′ sin θ cos ϕ+y′ sin θ sin ϕ+z′ cos θ.
Once the scattered fields are known, all quantities of interest regarding the light scattering can be calculated. For example, the time-averaged Poynting flow of the scattered field, or the scattered light intensity, at a far point (r,θ,ϕ), is given by:
Implementation and Validation of the Numerical Code
The computer code is written in FORTRAN 90 and the reported computations were carried out on Silicon Graphics Inc. (SGI) parallel computers. In most of the cases considered, the grid spacing used was λ/20 and the time step was equal to 0.8Δtc. Code validation was accomplished by comparison with Mie theory results. Unless noted otherwise, all calculations in this work use the following parameters: the plane wave is assumed incident along the positive z direction (i.e., θ=0 deg), the incident electric field is parallel to the x−z plane (i.e., ϕ=0 deg) and its amplitude is 1 V/m, the free-space wavelength is 0.6328 μm, the surrounding medium has a refractive index of 1.335, the cell has a refractive index of 1.591, and the position of the TF/SF interface is located 19 grid cells away from the boundary. The NTF interface is located 16 grid cells away from the boundary.
For the case of a plane wave incident on a uniform sphere, the scattering pattern can be described by the Mie theory. We have compared the FDTD results with those of the Mie theory for cases where the cell size and the scattering angle are varied. Figure 3 displays the scattering intensity versus the cell size at θ=0 deg. Figure 4 displays the scattering intensity versus the scattering angle. These two figures demonstrate the good agreement between our FDTD calculations with Mie theory results in both cell size and scattered angle scaling. Our calculations are valid over all scattering angles. To get the agreement at large angles approaching 180 deg, the effect of the numerical dispersion in the calculation of the far fields must be minimized, and the scattered field must be calculated to a high precision. One of the factors affecting the precision of the scattered field is multiscattereing, which refers to light that is scattered many times within the cell. A sufficiently long simulation time is necessary to ensure that multiscattering effects are properly taken into account.
The TF/SF approach implies that if there is no cell present, then in the TF domain the electromagnetic fields are purely incident fields, and in the SF domain the electromagnetic fields are equal to zero. Figure 5 displays the calculated value of Ex in x−z plane, with y being 2.133 μm from the FDTD boundary at time step nt=725. It is clearly shown that in the TF region (the center area of the plot), Ex has a sine wave form along the z axis and is constant along the x axis. In the SF region (border area of the plot), the maximum value of Ex is very small, and is close to 0. In this calculation, the TF/SF interface is located 0.474 μm from the FDTD boundary. To give a quantitative result, Fig. 6 displays the Ex waveform along the z axis, with x being 2.133 μm from the FDTD boundary. Other parameters in Fig. 6 are the same as those in Fig. 5. It is shown that in the SF region (two ends of the curve), Ex is of the order of 10−16 V/m, as is expected. The ramifications of the high accurary of the TF/SF implementation is the ability to reproduce light scattering at large angles where the intensity level is very small.
Features of Light Scattering from Cells
Over the past decade, a significant theoretical, computational, and experimental effort has been devoted to the exploration of the dependence of the scattering pattern on cell properties.3 4 12 13 14 In general, the cell size dominates the forward scattered intensity if the cell is surrounded by a medium with different refractive index. The nucleus is responsible for the scattering at slightly larger angles. Small organelles, or subcomponents of the cell, are responsible for the scattering that occurs at larger angles. Light scattering in the near backward direction is mostly due to the larger structures within the cell, such as the nucleus.
In cell diagnostics, if the size of cells to be diagnosed are significantly different, they are differentiable by only measuring the forward scattered (FS) light intensity.15 When the cells to be diagnosed are similar in size, such as monocytes and granulocytes in blood, a second parameter is required to distinguish them. The popular technique in present cytometry devices is to measure the side scattered (SS) light intensity, because the side scattered light intensity is very sensitive to the fine cell structure, especially organelles. We call this two-parameter (forward scattering and side scattering) cytometry. The cells are distinguished by displaying the measured data in 2-D coordinates of I SS and I FS , where subscripts FS and SS represent the forward scattered and side scattered light intensity, respectively. The different types of cells have different values of I SS and I FS , therefore they are in different positions in the I SS -I FS plot. If I SS and I FS for different types of cells are close to each other, these cells are difficult to distinguish. A good situation for cell differentiation is that for the same type of cells they are close to each other in the I SS -I FS plot, and for different types of cells they are far enough away from each other. Because the intensity of side scattered light and the differences in the side scattered intensity for different cells are all sensitive to the angle θ, the following question emerges: how to select an optimal angle for measuring the side scattered light so that the cells are most distinguishable? Due to the complexity of solving Maxwell’s equations, previous theories fail to give a clear answer to this question. Fortunately, our 3-D simulation is able to address this issue.
In two-parameter cytometry, the usual technique that is used to distinguish cells is to put the measured data on a x-y plot: the x axis represents the FS intensity and the y axis the SS intensity. The cells become distinguishable when these scattered data appear separated from one another in the x-y plot. We introduce a quantity δI to express the normalized distance between two cells in the x-y plot:δI value reflects the distinguishability of the cells. The larger the difference in δI between different types of cells, the better is the distinguishability of the cells. Because the dependence of the scattered intensity on θ is sensitive for large cells, the value of δI is also sensitive to θ. In practice, the measured light is taken over a range of angles in space, so we introduce δI, the area-averaged δI: I¯=∫ Ids/∫ ds, and the integration is over the cross section of the measured scattered light. In this work, the area to be averaged is indicated by a dimension δθ×δϕ. Later, we show that the selection of the size of the averaged area as defined by δθ and δϕ is also helpful for cell differentiation.
In the remainder of this section, we present examples that demonstrate how the features of light scattered from a cell can be determined, and how the strategy introduced earlier for determining the optimal angle for measuring side scattered light can be used in cell diagnostics based on two-parameter cytometry.
Case 1: Cells with the Same Structure but Different Refractive Index
For simplicity, we first assume that there are two types of cells with the same size and structure but with a different refractive index. Cells A and B are assumed to be spheres with a radius of 3.0 μm. The nucleus with a radius of 1.0 μm is located in the center of the cell. There are six subcomponents, each with a radius of 0.8 μm located on the coordinate axis 2.0 μm away from the original. The refractive index for cell A is 1.38 for nucleus and subcomponents and 1.36 for cytoplasm. The refractive index for cell B is 1.41 for the nucleus and subcomponents, and 1.375 for cytoplasm. The refractive index of the surrounding medium is assumed to be 1.33. In practice, the angle used in measuring side scattered light is in the vicinity of 90 deg or less. Our calculations are carried for θ less than 130 deg. Figure 7 displays the scattered intensity versus θ at ϕ=0 deg (parallel polarization, i.e., the scattered light is measured in the plane parallel to the incident electric field). Figure 8 displays the scattered intensity versus θ at ϕ=90 deg (perpendicular polarization, i.e., the scattered light is measured in the plane perpendicular to the incident electric field). It is shown that the scattered light patterns for two cells are similar and are barely distinguishable for most θ angles, except for a small range of angles around 90 deg for incident light with perpendicular polarization. However, the intensity value of fringer peaks of cell B is always larger than that of cell A. For angles other than where peaks occur, the value of the scattered light intensity of the two cells is very close. Therefore, the differentiation of these two cells in the x-y plot is strongly dependent on the angle where the scattered light is measured if the measured intensity is taken from a point, or from a very small area. These features can be better displayed by using the averaged δI parameter as shown in Figs. 9 and 10 for ϕ=0 and ϕ=90 deg, respectively. It is shown that when considering the dependence of δI on θ, δI is less sensitive to θ when the averaged angle increases. Therefore, an optimal δI is better suited for diagnostics use. These two figures show that these two types of cells are distinguishable over a large range of θ for both ϕ=0 deg and ϕ=90 deg. If ϕ=0 deg, the optimal angles lie in the small range of around 70 or 100 deg. If ϕ=90 deg, the optimal angles lie around 100 deg. In short, only cells with different refractive indices are easily distinguishable by using two-parameter cytometry at most angles. However, some optimal θ angles still exist for the case of polarized light. The optimal angle can be determined by our 3-D numerical simulation.
Case 2: Cells with the Same Size but Different Structure
In practice, there are often several types of cells in a sample, which are similar in size, e.g., lymphocyte, monocyte, and granulocyte in a blood sample. Here we assume that there are three types of cells with a radius of 4.5 μm. Cell A has a small 1-μm nucleus. Cell B has a large nucleus of 2.08 μm. Cell C has 72 organelles distributed within the cell. The radius of each organelle is assumed to be 0.5 μm. The refractive index of the nucleus and organelles is assumed to be 1.38. The refractive index of the cell cytoplasm is 1.36 and the refractive index of the surrounding medium is 1.33. The volume-averaged refractive index of cells B and C are the same. Figure 11 displays the scattered pattern from these three cells. It is shown that the patterns of scattered light are very similar and are barely distinguishable. Figures 12 and 13 display δI versus θ for parallel polarization and perpendicular polarization, respectively. It is shown that polarization properties can influence cell identification. In this case, when the side scattered light is observed at ϕ=0 deg, these three types of cells are only slightly distinguishable in two very narrow domains of θ around 70 and 100 deg, respectively. However, when the side scattered light is observed at ϕ=90 deg, the light scattered from these three type of cells is very distinguishable for θ in the range of 70 to 110 deg. Therefore, for polarized incident light, the selection of the angle ϕ is sometimes important. It can be concluded that cells with similar sizes but different structure are distinguishable only in some particular angular ranges. Furthermore, the polarization property of light is important for determining these optimal angles.
Case 3: Nonspherical Cells
In most cases, cells are modeled as spheres. This type of approximation can sometimes make analytical analysis possible. One such example is the Mie theory. However, many realistic cells are nonspherical. For nonspherical cells, the cell orientation will affect the pattern of scattered light. In other words, for the same type of cell, the scattered light pattern may be different due to the orientation of these cells. This will cause difficulty in cell differentiation. Minimizing the difference in the scattered light pattern for the same type of cells will be helpful for cell differentiation. Figure 14 displays the intensity patterns for an ellipsoidal cell when it is placed in one of two ways: the major axis along the z axis or transverse to the z axis. The refractive indices for cytoplasm and nucleus are 1.36 and 1.38, respectively. The radii of major and minor axes for the cell are 3.0 and 2.4 μm. The radii of major and minor axes for the nucleus are 1.0 and 0.8 μm. It is obvious that the scattered light patterns in these two situations are different, and for some particular angles the difference is very large. Figure 15 shows the distinguishability of these two positions: only in a few narrow angular ranges are the plots of these two positions very close, such as from 80 to 120 deg. Because one hopes that the scattering from the same type of cells is as similar as possible, the optimal angle for minimizing orientation effects can be selected based on the results in Fig. 15. It is necessary to understand that orientation of nonspherical cells also affects the cell distinguishability, and that the optimal angles for the minimizing orientation effect of nonspherical cells in the x-y plot can be found by our 3-D calculation.
A 3-D FDTD code is developed to simulate light scattering from complex cells and tissues in realistic situations. The code is also capable of accurately calculating backward scattered light. The simulation results provide an accurate description of the scattered light given the shape, structure, and electrical properties of cells. By introducing a normalized, area-averaged quantity δI, we present a method that can be used to find optimal angles where the cells are most distinguishable when using two-parameter cytometry. Our calculations also show that by selecting suitable scattering angles, cells with slight differences in their intrastructure can be distinguished with current cytometry techniques.
A detailed analysis of the angular dependence of the scattered light in two-parameter cytometry is carried out. It is shown that when the incident light is polarized, cell distinguishability is dependent on the incident light polarization. Some cell types are more easily distinguished by using parallel-polarized light, whereas other cell types are more easily distinguishable by using perpendicularly polarized light. This feature of the scattered light suggests that it is also important to select a suitable ϕ angle in cell diagnostics. In our calculations we have only displayed two cases: ϕ=0 deg and ϕ=90 deg. In practice, we can scan θ and ϕ over an entire angular space to find a set of optimal values for θ and ϕ. Cells with differences only in the refractive index are distinguishable at almost any angle; however, some optimal angles exist. Cells that have a similar size and an averaged refractive index, but differ in their internal structure, are distinguishable only over some special angles, and light polarization is important for cell distinguishability.
Due to the irregularity of cell shapes, cell orientation can also considerably affect the scattered light pattern. In cell diagnostics, the best angle for detecting the FS and SS light is the one in which the dependence on orientation is small for a given cell but large for other cell types. This can be deduced by careful examination of the δI plot, as shown in the last section.
Several examples of cases where a mixture of different cell types are present are also investigated. It is found that for all these cases, they are distinguishable by using two-parameter cytometry if suitable angles are selected to receive side scattered light. It is also found that for different cell mixtures, optimal angles for observing the side scattered light may be different. Fortunately, our 3-D simulations can give a very accurate prediction of the optimal choice for such angles. Examples given in this work only have the value of δI calculated in the vicinity of ϕ=0 deg and ϕ=90 deg. In practice, new cytometry devices are capable of measuring the side scattered light for ϕ varying over the range from 0 to 180 deg. Real cells typically have irregular shapes and their scattered light patterns might also be very different in the ϕ direction. A practical way to determine the optimal angle for distinguishing real cells is to calculate the value of δI over the entire measurable angular range and then find the angles where δI is maximum. Therefore, the combination of cytometry and 3-D simulations can be used to enhance the performance of current cytometry and to extend their capability to more subtle situations.
At present, the difficulty in cell diagnostics employing scattered light is two-fold: one is the difficulty in obtaining an accurate solution of Maxwell’s equations, in particular, in linking the measured scattered light to a particular cell; the other is the lack of practical knowledge of cellular structure and cellular electrical parameters (such as the refractive index and electric conductivity). The 3-D simulations provide a bridge that can connect the scattered light with the cell morphological and biophysical structure in realistic situations. Therefore, 3-D simulations can serve as a powerful tool not only in the accurate prediction of scattered light features, but also in the exploration and assessment of cell structure and its physical parameters.
The authors would like to thank M. Fujinaga of the MACI SGI Facility at the University of Alberta for his technical support. The work was supported by a strategic project of the Natural Sciences and Engineering Research Council of Canada and WestGrid computing resources. We would like to thank C. Backhouse and K. Singh for useful discussions.