Significant work has been done in characterizing the bulk optical absorption, scattering, and refractive properties of biological tissues separated into reasonable classes, such as skin, brain, bone, blood, etc.1 This work provides the foundational data for noninvasive optical analysis, which can, in theory, estimate the physiology and chemistry of an illuminated sample to a given degree of precision. That is, with the optical properties of the major constituents understood, it is possible to predict how light will interact with an anatomical part, such as a human finger, and extract from this interaction information about the blood and tissue. Of course, as has been summarized elsewhere,1 this is a difficult challenge since any optical property will change from person to person, from site to site on one person, and from time to time on any given site. It is, therefore, essential that detailed time- and position-sensitive data be modeled to understand what information can and cannot be extracted if noninvasive techniques are to be extended to measuring biological chemistry and/or physiology at increasingly finer detail.
One approach to understand how light propagates through tissues is to solve the diffusion equation as Steinke and Shepherd2 did with whole blood. However, the complex geometries that exist in tissue make this approach impractical without simplifying approximations. It has, therefore, been a common practice to use less constrained Monte Carlo simulations to model photon interactions in specific applications. We have adopted this approach, as have others,220.127.116.11.18.104.22.168.–12 because of its relative simplicity and accuracy.
Numerous studies involving Monte Carlo models have addressed the various geometries that are present within the human body, particularly with blood and tissue. Wang et al. created a multilayer grid model to explore diffuse reflectance and transmittance initiated by an infinitely narrow beam incident on the upper layer.8 Others have created more detailed models with more layers to explore skin reflectance1213.–14 and blood flow with layers.10,11 Some have explored the effects of complex boundaries.13 Still, others have modeled blood vessels within tissue,15 including using cylindrical models16 and by grouping the blood cells into discrete regions.17 This is a rich field, and we touch upon only a few of the many studies done within it.
As an exploration of the fundamental aspects of biological modeling itself, we have developed a Monte Carlo model to examine a homogeneous mix of basic biological material to understand better the complexities that may arise as structures are modeled in increasingly finer detail. We chose the finger for our own research interests, but the concepts and results presented in this article are general. We are interested in how near-infrared light transmits, particularly attenuates, through a semi-infinite slab of such a material. Others, such as Schmitt and Kumar,18 have used mathematical relationships and other physical principles applied to a mixed model to study how the scattering and back-scattering coefficients as well as the anisotropic factor vary based on the scattering contributions of the individual scatterers within a tissue as a whole. In the study reported here, we use Monte Carlo simulations to simulate light interaction with a mixture of two whole materials. Our aim is to find what effect mixing materials has on the attenuation of light within such a model as opposed to finding the spectral dependence of the scattering coefficients.
Monte Carlo Model
Our model begins with launching individual photon packets normally down into a specified material, where they interact in accordance with the material’s scattering and absorbing coefficients, and , respectively. The propagation distance, , between interactions is given as follows in terms of the interaction coefficient :6,7
The scattering of the photon packet is determined using directional cosines that include an anisotropic factor, and intensity is removed from the photon packet according to the absorption and scattering coefficients, as was done by Prahl et al.6
We set up our model in a Cartesian coordinate system, where the source is at the origin, and photon packets are launched in a collimated beam in the direction, normal to the interface. A circular detector 1 mm in diameter is placed at a location , where . We measure the intensity (termed “weight”) of the photon packets that reach the detector region. As we increase the thickness, , of the material, the transmitted intensity should follow a modified form of Beer’s Law for highly scattering materials given as19 We nominally run 50 simulations of varying thickness and from them, we derive the attenuation constant, . We keep the detector size and emitter intensity the same in all simulations so that the contribution to the scaling factor from the size of the detector is the same for all simulations and thus has no effect on the derived value of .
Methods and Procedures
There are cases where layered models are preferred over mixed models and vice versa. Blood and plasma are best represented in the finger by a mixed model since they travel together through blood vessels and veins. For epidermis and blood, a layered model is more useful in many regions, whereas a mixed model may be more accurate in other places. Although epidermis is a bloodless tissue,20 some areas of the epidermis layer are penetrated by capillaries, deviating from a layered model.21 A mixed model, although not completely accurate, approximates this penetration and is one step closer toward adding greater complexity to our models.
We chose to examine mixtures of epidermis, blood, and plasma because of their prevalence in the situations of interest to us21,22 and for simplicity in understanding how light attenuates differently for layered and mixed models. Since water makes up of plasma,23 we approximate plasma with water.
We ran semi-infinite models (finite in , infinite in and ) of pure material, layered material, and homogeneous mixed material using Mathematica. Illumination was from a collimated beam of monochromatic light at . We ran each simulation starting at a thickness of 0.010 mm and repeated the simulation at new thicknesses increased in increments of 0.100 mm up to a maximum thickness of 5.00 mm. We then plotted the transmitted intensity versus thickness and fit it with an exponential curve. The derived attenuation constant was then compared to given in Eq. (3) for the respective tissues.
Semi-Infinite Pure Tissue Model
To verify our approach, we first modeled semi-infinite pure tissue structures of epidermis, water, and whole blood, separately. These simulations successfully modeled the change in detected intensity versus the thickness of the sample for a modified form of Beer’s law in materials that exhibit high scattering properties [Eq. (3) and Fig. 1]. The derived values agreed with published values from diffusion models after they are normalized to our approach by dividing out the factor in Eq. (4) (see Table 1).
Absorption (μa), scattering (μs), transmission (μt), and bulk attenuation (α) coefficients and anisotropy factor (g) of epidermis (Ref. 19), water (Ref. 24), and whole blood [Hct=0.45−0.46, oxygenation >98% (Ref. 19)] at a wavelength of 810 nm. The values for αnormalized are the values calculated using Eq. (4), normalized by dividing out the 3. The values for αderived are the attenuation coefficients we obtained from fitting an exponential decay to the results from the simulations (see Fig. 1). Comparing αnormalized to αderived, we see that the two values agree quite well.
|Material||μa (mm−1)||μs (mm−1)||μt (mm−1)||g||αnormalized (mm−1)||αderived (mm−1)|
Semi-Infinite Layered Tissue Model
We next modeled layers of two different types of materials to compare against previously published results. This was accomplished by dividing up a semi-infinite volume into two separate layers, each composed of a different material. The interface between the two layers was given a specific coordinate, say , in the -plane (Fig. 2).
We modeled water and whole blood with water being the top layer and blood being the lower layer. We varied the percentage of each component from to and generated the same intensity versus thickness curves as earlier. We derived a value for as earlier and compared them with theory (Fig. 3).
One would expect there to be a linear relationship between the total attenuation and the volume fraction of each layer since the individual values mathematically add together. That is, if is the thickness of the layer of water, is the thickness of the layer of blood, is the attenuation coefficient of water, is the attenuation coefficient of blood, is the fraction of the total volume that is water, and is the fraction of the total volume that is blood, then multiplication of the two curves yieldsFigure 3 shows that the expected linear relationship between and the volume fraction of water is in excellent agreement with theory. We simulated a layered model of blood and epidermis as well and produced the expected linear relationship between volume fraction of either substance and the total attenuation coefficient in this case, too.
Semi-Infinite Homogeneous Mixed Tissue Model
We next modeled the light propagation through a homogeneous mix of two materials distributed randomly throughout the volume and packed loosely enough so that there was no overlap of interaction cross-sections. For both blood and water, we can accurately make this assumption. For epidermis, a tissue of tightly packed cells, this is not physiologically correct.20 Even so, we still consider mixed epidermis and blood. Those simulations can be thought of as gradually displacing the tightly packed epidermis cells with an increasing volume fraction of blood, thus reducing the packing factor that Schmitt and Kumar use in defining the of constituents of a particular tissue:182)]. For now, we neglect this packing factor in the epidermis and blood mixed model since we are simply interested in the effects on attenuation that a mixed model exhibits.
In deciding which material of a mixed medium was likely to interact with the photon packet, we formed a new , which was adjusted to account for the different interaction cross-sections of each type of material. Since each material has different scattering and absorption coefficients that are independent of volume, considering volume fraction alone is inaccurate. To illustrate why, compare the Monte Carlo photons to a bouncy ball moving around in a volume of bowling balls and ping-pong balls with enough space between the balls for the bouncy ball to move freely. The bowling balls represent particles of the material with higher scattering and absorption properties. The ping-pong balls represent particles of the material with lower scattering and absorption properties. Intuitively, it makes sense that as the bouncy ball moves around the volume, not only are the bowling balls going to absorb more energy than the ping-pong balls but it is also more probable that the bouncy ball will hit a bowling ball rather than a ping-pong ball due to the difference in sizes, or “interaction volumes,” of the two balls. (In this analogy, we relate a larger mass directly to a larger attenuation of the moving ball’s motion. In reality, a larger mass does not directly relate to attenuating the ball’s movement more, but we keep the analogy as an illustrative tool.)
Relating the simulation to the bowling-ball–ping-pong-ball analogy, if material 2 has higher scattering and absorption properties (e.g., the bowling balls) than material 1, then the photon packet is more likely to interact with material 2 than material 1. Thus, we rescaled the selection windows of each material to represent the fraction that each material takes of the whole interaction volume instead of the physical volume.
The rescaling goes as follows: let and have the same definitions as mentioned above, and let and be the interaction coefficients of materials 1 and 2, respectively. To rescale the volume fraction to the interaction volume of each material in the mixed model, the “cut-off” for the window in selecting between the two materials is
The propagation distance between interactions, , given by Eq. (1), is also affected by the difference in interaction cross-sections. In Eq. (1), the interaction coefficient, , corresponds to the overall interaction cross-section of a material. A larger corresponds to a larger interaction cross-section, and, in turn, a larger interaction volume. Similar to how Schmitt and Kumar rescaled the total scattering coefficient for the total scattering cross-section,18 and how Verkruysse et al. rescaled the absorption and scattering coefficients,17 we redefined as a weighted average of the interaction strengths of the two materials. Thus, the new , , has the form1).
Results and Discussion
After making the adjustments above, we ran the simulation for a mixture of water and blood and blood and epidermis. The results are in Figs. 4 and 5 and are plotted as versus volume fraction. These figures show that attenuation in a mixed medium does not have the linear relation with material fraction that is found in layered models (see Table 2).
The linear relation for in Beer’s law [Eq. (6)] is a consequence of the materials being physically and statistically independent. In the layered model, the materials are separated into distinct sections, where the propagation distance between interactions depends solely on the material the photon packet is in at any given moment [Eq. (1)]. Statistically, we can say that the probabilities of interaction of the different materials are independent of one another, as seen in the modified form of Beer’s law [Eq. (6)]. However, as soon as the materials are mixed together, the propagation distance between interactions is now dependent on all of the materials within the total interaction volume [Eq. (9)]. Therefore, the mixed model no longer yields statistical independence between the differing materials, and the linear relationship in no longer holds as it did in the layered model.
Comparison of total attenuation coefficient of a semi-infinite slab of 50% blood and 50% epidermis by volume derived from three different models: Monte Carlo multilayered (MCML), diffusion theory homogeneous mixed, and Monte Carlo homogeneous mixed.
|Homogeneous mix, diffusion theory (Ref. 17)||3.69|
|MC homogeneous mix||4.16|
In terms of the bowling-ball–ping-pong-ball analogy discussed in Sec. 3.3, a material with the larger , or larger interaction cross-section, would be represented by the bowling ball, whereas the one with the lower , or smaller interaction cross-section, would be represented by the ping-pong ball. If we took a mix of balls that was 50% bowling balls and 50% ping-pong balls by number and then separated them into two volumes, such as in the layered model, we would find that the bowling balls would occupy more volume than the ping-pong balls. Thus, would not prorate as 50% and 50% but rather have a higher weighting of than . As a result, the of the material with the higher , in this case, the bowling ball, would dominate the curve for , an argument that is also mentioned, but not fully characterized, by Schmitt and Kumar.18 The model with a mixture of water and blood is comparable to this analogy (Fig. 4). In this case, if blood is subscript “1” and water is subscript “2,” and .
In another situation, we could have while . Since of a given material is dependent on the anisotropic factor, , and is independent of , a situtation like this could exist if the anisotropic factor of material 2 is lower than that of material 1 [Eqs. (4) and (5)]. This scenario can be represented by large, light-weight plastic balls (material 1) and heavy shot-put balls (material 2). Since the plastic balls are larger, they have a higher interaction cross-section (), meaning it is much more likely for a photon packet to interact with them. On the other hand, since they are lightweight, they do not attenuate the intensity of the photon packet as much as the heavier shot-put ball would, which corresponds to a lower attenuation coefficient (). Since , will dominate the curve for . The model with a mixture of blood and epidermis is comparable to this analogy (Fig. 5).
The direction in which the curve for bends is determined by comparing of the two materials. The material with the higher will have a greater weighting in the total attenuation coefficient, . Thus, if the material with the greater also has a greater , then the curve for will bend above the linear relationship (Fig. 4). If the material with the greater has a smaller , then the curve will bend below the linear relationship (Fig. 5). How much the curve for deviates from the linear relationship depends on how much the values of for the two materials differ. In the case of the water–blood mix, (see Table 1), thus, the curve for bends far from the linear relationship (Fig. 4). On the other hand, in the blood–epidermis mixed model, , but they are of the same order of magnitude (Table 1). As a result, although the curve is more weighted to , the deviation from the linear relationship is small.
The approach of Verkruysse et al. for their analysis of the attenuation by blood in a port wine stain is similar to ours although there are some important differences to note.17 They use absorption and scattering coefficients and anisotropic factors that are linearly scaled according to the volume fraction of each material in the total mixture. The nonlinear aspects of increasing attenuation from greater vesicular blood volume is reflected in their correction factor . They then use these new mix coefficients to directly calculate the attenuation coefficient of a mixed model according to Eq. (4) (see Table 2). We linearly scaled the scattering and absorption coefficients according to the volume fraction of each material in the total mixture, which results in the transmission coefficient being linearly scaled according to volume fraction. However, we only use the rescaled transmission coefficient for the propagation distance and the selection window for which material to interact with since both of these parameters are dependent on all materials within the entire model. We did not rescale either or in the absorption factor (), and we do not rescale the anisotropic factor, . Our model is based on discrete events where a photon packet interacts with one material at each event. Both the absorption and anisotropic factors depend solely on the given event and, therefore, can only be based on the coefficients and factors of that one material.
Returning to the issue of the packing factor and its effects on introduced in Sec. 3.3, depending on whether a pure material is tightly packed or not, the nonlinear results could be mitigated or intensified. The results shown in Fig. 4 from the blood and water mixed model should not be affected as neither of the materials are tightly packed. However, in the case of the epidermis and blood mixed model, there should be some variations in the results once the packing factor is taken into account because tissue, such as epidermis, is considered to be tightly packed. By adding blood to a model of pure epidermis, we are essentially making the epidermis cells less tightly packed until we reach the model of pure blood. Reducing the packing factor, , of epidermis lowers the of epidermis, and in turn, lowers the of epidermis as well [Eqs. (2) and (7)]. Since of epidermis is already lower than that of blood (Table 1), the difference between the of the two materials becomes even greater, and thus, we should see the nonlinearities become more pronounced based on previous arguments.
We have developed a Monte Carlo homogeneous mixed model and have studied how light should attenuate in such a material. We have found that the linear scaling of attenuation coefficients found in layered models does not hold when materials are homogeneously mixed together. In cases, such as those involving blood or plasma, where there is no clear boundary between constituents or the constituents are mixed together in a fluid, it may be necessary to use a mixed model to derive the total attenuation coefficient with the required accuracy. This is especially true if the values of the materials differ by an order of magnitude or more; otherwise, a layered model may be sufficient.