Parallel operation of Monte Carlo simulations on a diverse network of computers

Accelerating fDOT image reconstruction based on path-history fluorescence Monte Carlo model by using three-level parallel architecture

The excessive time required by fluorescence diffuse optical tomography (fDOT) image reconstruction based on path-history fluorescence Monte Carlo model is its primary limiting factor. Herein, we present a method that accelerates fDOT image reconstruction. We employ three-level parallel architecture including multiple nodes in cluster, multiple cores in central processing unit (CPU), and multiple streaming multiprocessors in graphics processing unit (GPU). Different GPU memories are selectively used, the data-writing time is effectively eliminated, and the data transport per iteration is minimized. Simulation experiments demonstrated that this method can utilize general-purpose computing platforms to efficiently implement and accelerate fDOT image reconstruction, thus providing a practical means of using path-history-based fluorescence Monte Carlo model for fDOT imaging.

Mesh Optimization for Monte Carlo-Based Optical Tomography

Mesh-based Monte Carlo techniques for optical imaging allow for accurate modeling of light propagation in complex biological tissues. Recently, they have been developed within an efficient computational framework to be used as a forward model in optical tomography. However, commonly employed adaptive mesh discretization techniques have not yet been implemented for Monte Carlo based tomography. Herein, we propose a methodology to optimize the mesh discretization and analytically rescale the associated Jacobian based on the characteristics of the forward model. We demonstrate that this method maintains the accuracy of the forward model even in the case of temporal data sets while allowing for significant coarsening or refinement of the mesh.

Mesh-based Monte Carlo method using fast ray-tracing in Plücker coordinates

We describe a fast mesh-based Monte Carlo (MC) photon migration algorithm for static and time-resolved imaging in 3D complex media. Compared with previous works using voxel-based media discretization, a mesh-based approach can be more accurate in modeling targets with curved boundaries or locally refined structures. We implement an efficient ray-tracing technique using Plücker Coordinates. The Barycentric coordinates computed from Plücker-formed ray-tracing enables us to use linear Lagrange basis functions to model both media properties and fluence distribution, leading to further improvement in accuracy. The Plücker-coordinate ray-polygon intersection test can be extended to hexahedral or high-order elements. Excellent agreement is found when comparing mesh-based MC with the analytical diffusion model and 3D voxel-based MC code in both homogeneous and heterogeneous cases. Realistic time-resolved imaging results are observed for a complex human brain anatomy using mesh-based MC. We also include multi-threading support in the software and will port it to a graphics processing unit platform in the near future.

Full-field laser-Doppler imaging and its physiological significance for tissue blood perfusion

Using Monte Carlo simulations for a semi-infinite medium representing a skeletal muscle tissue, it is demonstrated that the zero- and first-order moments of the power spectrum for a representative pixel of a full-field laser-Doppler imager behave differently from classical laser-Doppler flowmetry. In particular, the zero-order moment has a very low sensitivity to tissue blood volume changes, and it becomes completely insensitive if the probability for a photon to interact with a moving red blood cell is above 0.05. It is shown that the loss in sensitivity is due to the strong forward scatter of the propagating photons in biological tissues (i.e., anisotropy factor g = 0.9). The first-order moment is linearly related to the root mean square of the red blood cell velocity (the Brownian component), and there is also a positive relationship with tissue blood volume. The most common physiological interpretation of the first-order moment is as tissue blood volume times expectation of the blood velocity (in probabilistic terms). In this sense, the use of the first-order moment appears to be a reasonable approach for qualitative real-time blood flow monitoring, but it does not allow us to obtain information on blood velocity or volume independently. Finally, it is shown that the spatial and temporal resolution trade-off imposed by the CMOS detectors, used in full-field laser-Doppler hardware, may lead to measurements that vary oppositely with the underlying physiological quantities. Further improvements on detectors' sampling rate will overcome this limitation.

Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration

General-purpose computing on graphics processing units (GPGPU) is shown to dramatically increase the speed of Monte Carlo simulations of photon migration. In a standard simulation of time-resolved photon migration in a semi-infinite geometry, the proposed methodology executed on a low-cost graphics processing unit (GPU) is a factor 1000 faster than simulation performed on a single standard processor. In addition, we address important technical aspects of GPU-based simulations of photon migration. The technique is expected to become a standard method in Monte Carlo simulations of photon migration.

Detection limits of multi-spectral optical imaging under the skin surface

The present work shows that the optical/biological information contained in a typical spectral image mainly reflects the properties of a small (conic like) volume of tissue situated vertically under each individual pixel. The objects appearing on a spectral image reasonably reproduce the correct geometrical shape and size (like a non-deformed shadow) of underlying inclusions of pathological tissue. The information contained in a spectral image comes from a depth that does not exceed approximately 2-3 mm. The number of photons that visit a given tissue voxel situated at a depth larger than approximately 2 mm represents less than the 1% of the total number of photons reaching the corresponding detection pixel (forming the image). A pathological inclusion (e.g. a pool of blood or vascular tumor) situated at a depth of approximately 0.5 mm with a thickness of 0.5 mm produces an image intensity contrast of approximately 5% (for images taken at wavelengths in the 600-1000 nm range) when compared to the normal skin background. The same inclusion at a depth of 20 microm provides a contrast decreasing from 55 to 20% with respect to an increase in wavelength. The dermis/hypodermis interface behaves as a partial barrier for the photons, limiting their access to deeper skin regions. The image contrast depends on the depth and the type of chromophore contained in the inclusion. An increase in the concentration of a given molecule may produce different contrast, independently of the depth, depending on the characteristics of the skin layer where this change occurs.

Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization

Monte Carlo (MC) based simulations of photon transport in living tissues have become the "gold standard" technique in biomedical optics. Three-dimensional (3D) voxel-based images are the natural way to represent human (and animal) tissues. It is generally believed that the combination of 3D images and MC based algorithms allows one to produce the most realistic models of photon propagation. In the present work, it is shown that this approach may lead to large errors in the MC data due to the "roughness" of the geometrical boundaries generated by the presence of the voxels. In particular, the computed intensity of the light detected on the tissue surface of a simple cubic tissue phantom may display errors from -80% to 120%. It is also shown that these errors depend in a complex manner on optical and geometrical parameters such as the interoptode distance, scattering coefficient, refractive index, etc. and on the degree of voxelization ("roughness") of the boundaries. It is concluded that if one wants to perform reliable 3D Monte Carlo simulations on complex geometries, such as human brain, skin or trabecular bone, it is necessary to introduce boundary meshing techniques or other equivalent procedures in the MC code to eliminate the deleterious effect of voxelization.

The use of the Henyey-Greenstein phase function in Monte Carlo simulations in biomedical optics

Monte Carlo (MC) simulations are often at the heart of the testing procedure in biomedical optics. One of the critical points in MC simulations is to define the new photon direction after each scattering event. One of the most popular solutions is to use the Henyey-Greenstein phase function or some linear combinations of it. In this note, we demonstrate that randomly generating the angle defining the new direction of a photon after a collision, by means of the Henyey-Greenstein phase function, is not equivalent to generating the cosine of this angle, as is classically done. In practice, it is demonstrated that for a nearly isotropic medium (asymmetry parameter g approximately 0) this discrepancy is not large, however for an anisotropic medium as is typically found in vivo (e.g. g = 0.98) the two methods give completely different results.

Simulation of emission tomography using grid middleware for distributed computing

SimSET is Monte Carlo simulation software for emission tomography. This paper describes a simple but effective scheme for parallel execution of SimSET using NetSolve, a client-server system for distributed computation. NetSolve (version 1.4.1) is "grid middleware" which enables a user (the client) to run specific computations remotely and simultaneously on a grid of networked computers (the servers). Since the servers do not have to be identical machines, computation may take place in a heterogeneous environment. To take advantage of diversity in machines and their workloads, a client-side scheduler was implemented for the Monte Carlo simulation. The scheduler partitions the total decay events by taking into account the inherent compute-speeds and recent average workloads, i.e., the scheduler assigns more decay events to processors expected to give faster service and fewer decay events to those expected to give slower service. When compute-speeds and sustained workloads are taken into account, the speed-up is essentially linear in the number of equivalent "maximum-service" processors. One modification in the SimSET code (version 2.6.2.3) was made to ensure that the total number of decay events specified by the user is maintained in the distributed simulation. No other modifications in the standard SimSET code were made. Each processor runs complete SimSET code for its assignment of decay events, independently of others running simultaneously. Empirical results are reported for simulation of a clinical-quality lung perfusion study.

Analysis of short-pulse laser photon transport through tissues for optical tomography

We describe a method for analyzing short-pulse laser propagation through tissues for the detection of tumors and inhomogeneities in tissues with the goal of developing a time-resolved optical tomography system. Traditional methods for analyzing photon transport in tissues usually involve the parabolic or diffusion approximation, which implies infinite speed of propagation of the optical signal. To overcome such limitations we calculate the transmitted and reflected intensity distributions, using the damped-wave hyperbolic P(1) and the discrete-ordinates methods, for a wide range of laser, tissue, and tumor parameters. The results are compared with the parabolic diffusion P(1) approximation.

Relevance of accurate Monte Carlo modeling in nuclear medical imaging

Monte Carlo techniques have become popular in different areas of medical physics with advantage of powerful computing systems. In particular, they have been extensively applied to simulate processes involving random behavior and to quantify physical parameters that are difficult or even impossible to calculate by experimental measurements. Recent nuclear medical imaging innovations such as single-photon emission computed tomography (SPECT), positron emission tomography (PET), and multiple emission tomography (MET) are ideal for Monte Carlo modeling techniques because of the stochastic nature of radiation emission, transport and detection processes. Factors which have contributed to the wider use include improved models of radiation transport processes, the practicality of application with the development of acceleration schemes and the improved speed of computers. In this paper we present a derivation and methodological basis for this approach and critically review their areas of application in nuclear imaging. An overview of existing simulation programs is provided and illustrated with examples of some useful features of such sophisticated tools in connection with common computing facilities and more powerful multiple-processor parallel processing systems. Current and future trends in the field are also discussed.