PLATFORM FOR COMBINED MONTE CARLO AND K-WAVE SIMULATIONS OF OPTOACOUSTIC
IMAGES FOR BLOOD SATURATION MAPPING
ALEKSANDR KHILOV, VALERIYA PEREKATOVA, DARIA KURAKINA, AND MIKHAIL KIRILLIN
Institute of Applied Physics RAS, Nizhny Novgorod, Russia [email protected]
ABSTRACT
Optoacoustic (OA) imaging [1-5] is a modern non-invasive hybrid visualization technique based on tissue illumination with short laser pulses and further detection of acoustic waves generated by local heating of the tissue caused by the absorption of probing radiation. One of the most promising directions of OA imaging application is the mapping of blood oxygen saturation StO2 based on OA spectroscopy data [6, 7] owing to significant difference in oxy-and deoxyhemoglobin absorption spectra. Biological tissues are optically inhomogeneous, and their optical properties are a priori unknown, which requires development of customized approaches for extraction of physiological parameters. Current trend consists in application of machine learning techniques [8, 9] using large sets of synthetic OA data, which requires numerical solutions of both optical and acoustic problems.
We report on the platform for numerical simulation of OA data based on three-dimensional Monte Carlo (MC) simulation [10-12] of light propagation in biological tissue for the solution of optical problem combined with the calculation of acoustic wave propagation with k-Wave toolbox [13]. The developed platform was employed for calculations of reconstructed OA images for several probing wavelengths, 532, 658 and 1064 nm, which are common for spectral OA imaging aimed at StO2 reconstruction [14].
The full size of a tissue sample considered in simulations was 20x20x10 mm3. The light source was considered to be a plane wave illuminating the sample from the top. For each calculation, a total number of 107 photons were launched into the medium perpendicular to the tissue surface and used to calculate absorption maps for each of three probing wavelengths. Calculated absorption maps (photon weight per volume unit) were normalized to initial intensity (photon weight per area unit) thus resulting to maps dimension of reversed length (mm-1). The simulations were performed for simplified skin tissue model (Fig. 1a-c), containing 4 blood vessels (v1-v4) with different diameters (varied in range of 0.05-1 mm) and embedding depths (varied in range of 0.1-2 mm). One configuration contained 4 vessels with different parameters resulting in a total of 50 different configurations. The background tissue is considered to be flat with homogenous optical properties. Blood oxygen saturation in blood vessels was varied from 0.1 to 1 with the step of 0.1. The blood optical properties were calculated based on absorption and scattering spectra for oxygenated and deoxygenated whole blood taken from paper [15]. Optical properties of the background tissue were taken for human skin in vivo from paper [16]. Thus, the total number of calculated absorption maps for different 50 geometries, 3 probing wavelengths and 10 oxygenation values amounted 1500.
'A Ä A b) 1 mm
1 mm
h)
1 mm
~~ 1
y
1
0.8 0.6 0.4
0.2 x10! 8 6 4
2 0 -2
• C) 1 mm
1
f) 1 mm
I
0.3
0.1
150
Figure 1: Initial normalized pressure distributions (units of mm'1) for OA modeling obtained from MC (a,b,c), k-Wave toolbox simulated raw (d,e,f) and reconstructed (g,h,i) OA B-scans. All bars are 1 mm
Figures 1a-c show central cross-sections of three-dimensional maps of initial pressure distributions obtained by multiplication of absorbed light dose by Gruneisen parameter in blood vessels for their various diameters and embedding depths for probing wavelengths of 532 nm (Fig 1a), 658 nm (Fig. 1b) and 1064 nm (Fig. 1c). At the second stage, based on the obtained initial pressure distributions, the OA microscopy B-scans (Fig. 1d-f) were obtained in three-dimensional geometry by the solution of the forward acoustic problem using k-Wave [13] for each detector position in XZ-plane. The ultrasonic detector was supposed to be a spherically focused antenna with the parameters close to those for real detector employed in an OA microscope [17-19]. After k-Wave modeling, each OA microscopy B-scan was processed with the synthetic aperture focusing technique described in paper [20]. The reconstructed B-scans at three wavelengths are shown in Fig.1 g-i.
Figures 2a-c show the dependencies of ratio of OA signal intensities registered at different probing wavelengths and averaged over each of four blood vessels (v1-v4 in Fig. 1a) on blood oxygen saturation, which were obtained in
32
numerical simulation. The signal ratios show monotonous dependence on StO2 due to monotonous dependence on this parameter. Meanwhile, they show different trends due to difference in oxy- and deoxyhemoglobin absorption spectra. The demonstrated monotonous dependence of signal ratio (Fig.2) indicates the possibility of employing the ratiometric approach in StO2 reconstruction. The dependence of ratio of OA signals corresponding to probing wavelengths of 658 nm and 1064 nm on blood oxygen saturation demonstrates the smallest variance for all diameters and depths of vessels and can be employed for a simple estimation of StO2 value.
532nm/1064 nm
532rim/658 nm
658nm/1064 nm
Figure2: Dependencies of the ratios of calculated OA signal of blood vessel on blood oxygenation level @ 532 nm / @1064 nm (a), @
532 nm / @658 nm (b) and @ 658 nm / @1064 nm (c)
Reconstruction of blood oxygen saturation was based on inverse dependence of blood oxygen saturation vs ratio of intensities. Since the ratio of intensities corresponding to 658 nm and 1064 nm is less dependent on blood oxygen saturation, the closest trends (v2, v3, v4 in Fig. 2c) were averaged to construct inverse function StO2(p658/p1064). After that, calculated maps of 658/1064 ratio were employed for the reconstruction of blood oxygen saturation based on inverse function StO2(p658/p1064).
The results of StO2 reconstruction for several simple geometry cases are shown in Fig. 3. Each figure contains schematic of vessel location with real StO2 value (Fig. 3a,e,i), calculated OA signals ratio maps in XZ projections (Fig.3b,f,j), projections of full reconstructed maps of StO2 (Fig.3c,j,k) and reconstructed StO2 values in vessels only (Fig.3d,h,l).
Figure3: Vessel masks, calculated 658nm/1064 nm ratio, reconstructed blood oxygen saturation, and reconstructed blood oxygen saturation in vessels only (from left to right) for initial blood oxygen saturation 0.6 (top row), 0.8 (middle row) and 1 (bottom row). All
bars are 1 mm
As can be seen from Fig. 3, the proposed simple voxel-by-voxel reconstruction algorithm of blood oxygen saturation value based on the OA signals ratio provides with adequate estimation of blood oxygen saturation in vessels with no a priori knowledge of vessel embedding depth required. For true StO2 values >0.8 revealed precision <15%. Moreover, depth-dependent StO2 reconstruction is expected to provide even higher accuracy. In this connection, employing a more advanced, however, time-consuming approach to OA image simulations in combination with machine learning techniques instead of the considered simple algorithm has high potential in accurate mapping of blood oxygenation based on spectral OA data.
ACKNOWLEDGEMENTS
The study is supported by Russian Science Foundation (project #22-29-00074). REFERENCES
[1] V. Perekatova, M. Kirillin, P. Subochev, A. Kurnikov, A. Khilov, A. Orlova, D. Yuzhakova, I. Turchin, Quantification of microvasculature parameters based on optoacoustic angiography data, Laser Phys. Lett. 18, p. 035602, 2021.
[2] V. Perekatova, S. Nemirova, A. Orlova, M. Kirillin, A. Kurnikov, K. Pavlova, A. Khilov, A. Kovalchuk, P. Subochev, Three-dimensional dual-wavelength optoacoustic angiography reveals arteriovenous anastomoses, Laser Phys. Lett. 18, p. 045601, 2021.
[3] L.V. Wang, S. Hu, Photoacoustic tomography: in vivo imaging from organelles to organs, Science 335, p. 1458-62, 2012.
[4] M. Xu, L.V. Wang, Photoacoustic imaging in biomedicine, Rev. Sci. Instrum. 77, p. 041101, 2006.
[5] P. Beard, Biomedical photoacoustic imaging, Interface Focus 1, p. 602-31, 2011.
[6] V. Perekatova, P. Subochev, M. Kirillin, E. Sergeeva, D. Kurakina, A. Orlova, A. Postnikova, I. Turchin, Quantitative techniques for extraction of blood oxygenation from multispectral optoacoustic measurements Laser Phys. Lett. 16, p. 116201, 2019.
[7] B.T. Cox, J.G. Laufer, P.C. Beard, S.R. Arridge, Quantitative spectroscopic photoacoustic imaging: a review, J. Biomed. Opt. 17, p. 061202, 2012.
[8] J. Grohl, M. Schellenberg, K. Dreher, L. Maier-Hein, Deep learning for biomedical photoacoustic imaging: A review, Photoacoustics 22, p. 100241, 2021.
[9] C. Yang, H. Lan, F. Gao, F. Gao, Review of deep learning for photoacoustic imaging, Photoacoustics 21, p. 100215, 2021.
[10] M. Kirillin, V. Perekatova, I. Turchin, P. Subochev, Fluence compensation in raster-scan optoacoustic angiography, Photoacoustics 8, p. 59-67, 2017.
[11] M. Kirillin, A. Khilov, D. Kurakina, A. Orlova, V. Perekatova, V. Shishkova, A. Malygina, A. Mironycheva, I. Shlivko, S. Gamayunov, Dual-Wavelength Fluorescence Monitoring of Photodynamic Therapy: From Analytical Models to Clinical Studies, Cancers 13, p. 5807, 2021.
[12] A.V. Gorshkov, M.Yu. Kirillin, Acceleration of Monte Carlo simulation of photon migration in complex heterogeneous media using Intel many-integrated core architecture, J. Biomed. Opt. 20, p. 085002, 2015.
[13] B.E. Treeby, B.T. Cox, k-Wave: MATLAB toolbox for the simulation and reconstruction of photoacoustic wave fields, J. Biomed. Opt. 15, p. 021314, 2010.
[14] V. Perekatova, P. Subochev, M. Kleshnin, I. Turchin, Optimal wavelengths for optoacoustic measurements of blood oxygen saturation in biological tissues, Biomed. Opt. Express 7, p. 3979-95, 2016.
[15] N. Bosschaart, G.J. Edelman, M.C. Aalders, T.G. van Leeuwen, D.J. Faber, A literature review and novel theoretical approach on the optical properties of whole blood, Lasers Med. Sci. 29, p. 453-79, 2014.
[16] T. Kono, J. Yamada, In vivo measurement of optical properties of human skin for 450-800 nm and 950-1600 nm wavelengths, Int. J. Thermophys. 40, p. 1-14, 2019.
[17] I. Turchin, S. Bano, M. Kirillin, A. Orlova, V. Perekatova, V. Plekhanov, E. Sergeeva, D. Kurakina, A. Khilov, A. Kurnikov, P. Subochev, M. Shirmanova, A. Komarova, D. Yuzhakova, A. Gavrina, S. Mallidi, T. Hasan, Combined Fluorescence and Optoacoustic Imaging for Monitoring Treatments against CT26 Tumors with Photoactivatable Liposomes, Cancers 14, p. 197, 2021.
[18] V. Perekatova, M. Kirillin, S. Nemirova, A. Orlova, A. Kurnikov, A. Khilov, K. Pavlova, V. Kazakov, V. Vildanov, I. Turchin, P. Subochev, Quantitative Characterization of Age-Related Changes in Peripheral Vessels of a Human Palm Using Raster-Scan Optoacoustic Angiography, Photonics 9, p. 482, 2022.
[19] P. Subochev, Cost-effective imaging of optoacoustic pressure, ultrasonic scattering, and optical diffuse reflectance with improved resolution and speed, Opt. Lett. 41, p. 1006-9, 2016.
[20] V.V. Perekatova, M.Yu. Kirillin, I.V. Turchin, P.V. Subochev, Combination of virtual point detector concept and fluence compensation in acoustic resolution photoacoustic microscopy, J. Biomed. Optics 23, p. 091414, 2018.