Photoacoustic Tomography predictive model applied to basocellular carcinoma for several sensor geometries

A predictive model for Photoacoustic Tomography images is presented and applied to basocellular carcinoma. It implements optical propagation by a Monte Carlo approach, initial pressure distribution and image reconstruction by the k-space method. OCIS codes: (110.5120) Photoacoustic imaging; (170.1870) Dermatology; (170.3660) Light propagation in tissues


Introduction
Clinical diagnosis can be greatly improved in the next years by Photoacoustic Tomography [1].Photoacoustic Tomography is able to improve the penetration depth of optical techniques [2], while preserving high resolution images.This technique is non-invasive, as it employs a combination of an optical source with the internal generation of acoustic waves, which are externally detected.Photoacoustic Tomography (PAT) is based on the photoacoustic effect.The object, biological tissue in this case, is irradiated by a short pulsed laser beam under certain confinement conditions.Part of the light is absorbed by the object and partially converted into heat.The heat is then converted into a pressure rise by means of thermoelastic expansion.The pressure rise propagates as an ultrasonic wave, is detected by transducers and employed to form an image [1].Available optical ballistic imaging modalities, such as confocal microscopy, two-photon microscopy, or optical coherence tomography cannot penetrate into biological tissues beyond the soft depth limit, at about 1-2 mm [2].Diffuse optical tomography is able to go further than this limit, but at the expense of a reduced resolution.PAT contributes to obtain high-resolution images beyond this limit.Apart from this advantage, PAT has a high endogenous contrast given by local optical absorption.
Optical or photoacoustic techniques are particularly useful in dermatological diseases.The most common type of cancer in Western countries is nonmelanoma skin cancer, which accounts for approximately one third of all diagnosed cancers [3].It represents around 80% of the total amount of cutaneous cancers, with an apparition rate 20 times that of the melanoma.Among all these cases of nonmelanoma skin cancer, Basal Cell Carcinoma (BCC) represents more than one half.
When evaluating the performance and validity of Photoacoustic Tomography, a deep knowledge of the processes involved is required.In this context predictive tools are highly desired, as the efficacy for particular medical applications could be estimated [4].The limitations of the technique could also be established.
In this work a predictive model for photoacoustic images, based on a Monte Carlo method for optical propagation, and the k-space method for photoacoustic reconstruction, is developed.The model is applied to dermatological tumoral tissue, particularly infiltrative BCC.First in section 2 the predictive analysis model, including optical propagation, initial pressure distribution and image reconstruction by the k-space method, is described.Afterwards section 3 shows the application to infiltrative BCC diagnosis and the results obtained.

Photoacoustic Tomography predictive model
Photoacoustic Tomography analysis requires first an estimation of optical propagation in the biological tissue.Analytic solutions for realistic scenarios are usually complex, or even not affordable.Numerical techniques are much more recommended in this kind of applications that involve turbid media like biological tissues.Approximations based on Monte Carlo (MC) approaches are quite employed [2].The step size and angular scattering deviation at each event are obtained by statistically sampling the probability distributions.As the number of photons increases, the distribution of all paths produces an increasingly accurate approximation of the problem.The actual number of photons required for a realistic result depends on the details of the simulation.
In this work the Monte Carlo simulation package MCML is used [5,6].The first step is to define the optical properties [2] and thickness of each layer of the biological tissue.After that, the laser beam is defined, mainly spatial profile, radius and total energy, and the convolved data are obtained.
Once the optical propagation in the biological tissue is known, the photoacoustic effect must be taken into account.Two conditions must be met for the optimal generation of optoacoustic signals, thermal and stress confinement [7].Under these conditions the initial pressure increase after laser excitation can be expressed as follows: Where β(K −1 ) is the thermal coefficient of volumetric expansion, γ(Pa −1 ) is the isothermal compressibility coefficient, ρ � kg m 3 � is the tissue density, and C V � J Kg � denotes the specific heat capacity at constant temperature.
The Grüneisen parameter is relevant for the photoacoustic effect.It is dimensionless and indicates the fraction of incident optical energy converted into acoustic energy by the photoacoustic effect.It can be defined as [7]: And consequently, relating the initial pressure distribution with the optical propagation in the biological tissue: Where F(J/cm 2 ) is the fluence rate of the laser source.After the initial pressure calculation the acoustic wave starts to propagate.The analysis of the acoustic propagation is complex, and can be solved numerically.One approach consists of using the k-space method to solve the system of coupled acoustic equations [8].It allows the reduction of time steps required for accurate simulations, and as a consequence of the memory required.This approach combines the spectral calculation of spatial derivatives with a temporal propagator expressed in the spatial frequency domain or k-space.In a standard finite difference scheme, spatial gradients are computed locally based on the function values at neighbouring grid points.In the simplest case, the gradient of the field can be estimated using linear interpolation.A better estimate of the gradient can be obtained by fitting a higher-order polynomial to a greater number of grid points and calculating the derivative of the polynomial.The Fourier collocation spectral method takes this idea further and fits a Fourier series to all the data.There are two significant advantages of using Fourier series.First, the amplitudes of the Fourier components can be calculated efficiently using the Fast Fourier Transform (FFT).Second, the base functions are sinusoidal, so only two grid points (or nodes) per wavelength are theoretically required, rather than the six to ten required in other methods.
In this work the k-space solver will be based on k-wave [9], an acoustic toolbox for Matlab © for photoacoustic image reconstruction.Several sensor geometries will be considered, from the standard linear, L-shaped or iterative, to the newly implemented U-shaped or squared.

Application to basocellular carcinoma
The previous model was applied to carcinoma diagnosis, particularly to dermatological infiltrative BCC.The BCC was considered to be 1.3 mm deep, and to be over three layers of healthy skin.The optical properties of the biological tissue can be found elsewhere [10].The laser radius was considered to be 0.25 mm, with 200 mW/cm 2 irradiance at a wavelength of 520 nm.Optical propagation was calculated by means of the Monte Carlo method.Afterwards the initial pressure was computed, Fig. 1, and the acoustic propagation and image reconstruction was calculated by k-wave.Several sensor geometries were tested, mainly linear, L-shaped, iterative, U-shaped and squared.Fig. 2 shows the U-shaped geometry, and Fig. 3 the squared geometry.In order to test the quality of the reconstruction, Fig. 4 shows the real and estimated initial pressure for each sensor geometry.As it can be easily appreciated, especially in Fig. 4, both the iterated and the squared reconstruction approaches provide the best results for the estimation of initial pressure, what is directly related with carcinoma diagnosis.

Fig. 1 .
Fig. 1.Radial and depth initial pressure distribution in the tissue.Fig. 2. Reconstructed initial pressure distribution with an U-shaped sensor.

Fig. 3 .
Fig. 3. Reconstructed initial pressure distribution with a squared sensor.Fig. 4. Comparison between calculated and reconstructed initial pressure distributions for different sensor geometries.