CC BY
50
Kolesenkov Alexander Nikolaevich, candidate of technical sciences, docent, sk62@,mail.ru, Russia, Ryazan, Ryazan state radio engineering University,
Fetisov Dmitry Vadimovich, postgraduate, morzitko@,gmail. com, Russia, Ryazan, Ryazan state radio engineering University
УДК 004
МОНТЕ-КАРЛО МОДЕЛИРОВАНИЕ СТРУКТУРНЫХ ОКТ ИЗОБРАЖЕНИЙ КОЖИ ЧЕЛОВЕКА IN VIVO С ИСПОЛЬЗОВАНИЕМ ЭКСПЕРИМЕНТАЛЬНЫХ В-СКАНОВ И РАСПРЕДЕЛЕНИЯ ОПТИЧЕСКИХ СВОЙСТВ
С.Н. Абдулкарим, Д.А. Петров, С.Г. Проскурин
Описано моделирование структурных ОКТ (оптическая когерентная томография) изображений методом статистических испытаний с использованием экспериментальных изображений реальных биообъектов. Биологический объект рассматривается как совокупность элементов 3D, которые позволяют моделировать среду, структура которой не может быть описана аналитически. Каждый воксел характеризуется своим показателем преломления и анизотропии, рассеянием и поглощением. Сканирование внутренней структуры биообъекта используются для восстановления моделированного изображения вместо аналитического представления о геометрии границ неоднородностей. Эффективность описанного метода проверена путем сравнения моделируемых и экспериментально полученных А- и В-сканов.
Ключевые слова: оптическая когерентная томография, глубина когерентного зондирования, человеческая кожа, моделирование методом Монте-Карло, воксельная модель среды.
Оптическая когерентная томография (ОКТ) - неинвазивный метод медицинской визуализации основанный на принципах низкокогерентной интерферометрии. Благодаря высокой разрешающей способности и безопасности пациентов этот метод в настоящее время используется в различных областях медицины, таких как офтальмология, стоматология и дерматологии [1]. Однако ОКТ имеет ряд существенных недостатков, которые ограничивают его эффективность. Главным ограничивающим фактором является глубина когерентного зондирования (ГКЗ), которая из-за большого рассеяния оптического излучения, как правило, не превышает 1-2 мм в зависимости от исследуемого биообъекта [2-5]. Понимание механизмов
227
рассеяния света и упрощение процесса получения структурных изображений необходимо, чтобы улучшить контраст и увеличить ГКЗ. На данный момент существует лишь несколько статей, описывающих Монте-Карло моделирование и реконструкции ОКТ изображений. Все они, как правило, основаны на алгоритме статистических испытаний, который реализует пошаговый перенос внутри сильно рассеивающего исследуемого объекта. Длина шага определяется как:
s = (1)
m«+ m ^
где ma и m^- коэффициенты поглощения и рассеяния среды, соответственно, X - случайное число, равномерно распределены между 0 и 1.
Фазовая функции Хеньи-Гринштейна обычно используется для определения угла рассеяния фотона (урав. 1):
1 - g 2
P(cos 0) =-—-— (2)
2(1 + g2 - 2g cos 0) 2
где 0 - угол рассеяния, g - параметр анизотропии.
После этого, определяется статистический вес фотонного пакета, и весь процесс повторяется до тех пор когда фотон достигает границы объекта или его «вес» стал ниже установленного порога [6-13]. Из-за статистической природы метода Монте-Карло время вычисления ОКТ изображений становиться очень длинным (2). Реальные биологические объекты, как правило, имеют сложную внутреннюю структуру, которая не может быть описана с помощью простых аналитических формул или равномерного геометрического распределения оптических свойств.
В литературе описаны различные способы задания геометрии неод-нородностей и построения объекта. Самый простой - это способ описания биообъекта в виде ряда слоев с различными оптическими свойствами и прямых границ между ними [7, 8]. Этот подход предлагает быстрый и простой способ, чтобы определить точку пересечения луча и границы, где размер шага фотона, s^, можно рассчитать следующим образом:
Sb = ^, (3)
uz
где z - текущее положение фотона на оси z, zo - граница координат, а uz - направление движения фотона. После этого, определяют расстояние между рассеяниями по сравнению с исходным размером шага. После этого, определяют, продолжает ли фотон мигрировать в текущем слое или перейдет к следующему (3). Однако реальный биомедицинский объект может быть лишь приближенно рассматриваться, как ряд слоев с таким типом границ. Оптическая диагностика, как правило, имеет дело со сложной геометрией распределения оптических неоднородностей.
Другой способ задания границ заключается в использовании аналитических уравнений [9, 10]. Такой подход позволяет представить сложные структуры, такие как эпидермис и дерма слоёв кожи. В этом случае, вычисление расстояний зависит от сложности граничных условий, что в большинстве случаев означает, что необходимо использовать численное приближение, например Метод Ньютона, который делает моделирование ещё более медленным.
В некоторых случаях можно сделать моделирование более простым и быстрым, используя аналитический подход. Он может быть применён, когда исследуемые неоднородности внутри объекта задаются простой геометрии, например, параллелограммом или эллипсоидом [11]. Этот метод может быть использован в сочетании с аналитическим описанием, однако, в случае, когда внутри тканей существует множество геометрических объектов, время вычисления опять становится значительным.
Можно использовать более гибкий метод, в котором реализуется полигон сетки, представляющий распределение оптических свойств для моделирования методом Монте-Карло [12, 13]. В этом случае геометрия границ представлена как последовательность объемов тетраэдров, каждая из которых имеет свои оптические свойства. Расстояние до точки пересечения луча и плоскости тетраэдра вычисляется по формуле:
П.р + (
=--, (4)
пм
где П - нормальный единичный вектор в плоскости пересечения, м и р-положение и направление фотона.
Слабым местом этого метода является большое число ячеек, необходимых для представления объекта и, следовательно, большое количество облучения границ пересечений, которые требуют слежение за траекторией миграции фотонов в ткани. Кроме того, существуют и другие подходы, которые реализуют сегментацию моделируемой структуры биообъекта для метода вокселизации границ объектов (4). Они имеют те же проблемы с увеличением времени необходимых расчетов для учета пересечения границ рассеянным излучением.
Чтобы сделать моделирование одного А-скана для реконструкции ОКТ-изображений с использованием метода Монте-Карло требуется ис-
7 9
пользовать 10-109 пакетов фотонов. Поэтому время моделирования может быть значительно увеличено, особенно в случае, если количество получаемых А-сканов имеет важное значение. Поскольку наибольшее количество фотонов рассеивается внутри объекта и не достигает детектора, искусственная функция рассеяния может быть реализована для прямых пакетов фотонов в направлении зонда. В простейшем случае, после первых собы-
тий рассеяния фотона рассчитывается направление рассеяния, (см. (2)), чтобы определить новое направление и его «вес» [8]. Модернизация этого подхода была осуществлена с помощью функции определённой плотности вероятности для первого рассеяния, которое перенаправляет фотон в направлении положения детектора [13]. Для расчёта дальнейшей траектории движения, рассеяния фотонов определяется заранее с некоторой заданной вероятностью и также может направляться к детектору. Применение этого метода должно дать десятикратное увеличение скорости моделирования, однако, до сих пор проверка этого метода с использованием структурных ОКТ изображений не была реализована.
В этой статье мы описываем гибкий метод моделирования ОКТ изображений, который основан на вокселном подходе к моделируемому объекту, который реконструируется с использованием экспериментально полученных структурных изображений. Валидация методики проводится с использованием реальных ОКТ изображения поверхностных слоёв кожи человечка in vivo без применения метода снижения дисперсии.
Геометрия моделируемого объекта. Подход вокселизации к представлению объекта позволяет описать объект с произвольной и сложной геометрией. Разрешение распределения оптических свойств зависит от количества вокселей и их размеров (рис. 1). По мере того как число типов тканей, как правило, ограничены, каждый из вокселей присвоен номер индекса соответствующего набора оптических свойств.
Рис. 1. Вокселизация трехмерного объекта
Облучение - воксельная краевая задача пересечения решается с помощью алгоритма Смита:
dbr
z + dz - z с -z-- if uz > 0
uz
z-z
с
u
if uz < 0
dby =
У + dy - y-
u
if uy > 0
y
У - У- •
(5)
u
if uy < 0
y
db
x
X + dx Xс if ux > 0
xz X X
с
x.
if ux < 0
где dbx,dby,dbz - расстояния до воксельной границы; хс,ус,zс - текущее положение фотона; dx, dy, dz - воксельная длина ребра; z, у, х - координаты нижней границы воксела.
Каждые расстояния по сравнению с более низкими значениями представляют собой размер шага до точки пересечения.
Экспериментальные В-сканы используются для моделирования объекта с использованием соответствующего алгоритма обработки, который может быть реализован с помощью обнаружения слоя [14]. Неровное изображение, где каждый цвет представляет собой биологическую ткань может быть преобразована в 3D объект (^скан) для моделирования с использованием подхода вокселизации (рис. 2).
Рис. 2. Экспериментальный В-скан (2D структурное ОКТ изображение) кожи пальца человека in vivo преобразуется в С-скан
для дальнейшего моделирования
Моделирование структурных изображений. Основой представленного алгоритма моделирования структурных изображений является классическим методом МКПФ (Монте-Карло для переноса фотонов) с некоторыми изменениями, сделанными намеренно, чтобы реализовать реконструкцию 2Б структурного изображения (В-скана). Блок-схема алгоритма, представлена на рисунке 3.
Фотонный пакет запущен в исходной точке в соответствии с гауссовым распределением вероятностей формы пучка. Каждый раз, когда фотон мигрирует внутри объекта алгоритм Смита используется для проверки текущего состояния (урав. 5), размер шага определяется по (урав. 1). Для вычисления длины шага, вместо полного коэффициента экстинкции мы используем только коэффициент рассеяния. Поглощение определяется законом Ламберта-Бугера:
Рис. 3. Алгоритм структурного моделирования ОКТ изображения
232
Когда фотон выходит из ткани внутри зонда с радиусом й интенсивности каждого пиксела, содержащегося в текущем А-скана уменьшается в соответствии с:
1 'i,j=
■Jw exp
r2z - L
lc у
005
2p(2z - L) l
if uz >0; r < d
(6)
0, otherwise
где Ж - статистический вес фотонного пакета; 1с - длина когерентности; 2 г и Ь - справочная и образец длины оптического пути, соответственно; 0 - детектор угла приема; i - ширина А-скана индекс; у - глубина А-скана индекс; 1 - длина волны источника; г - расстояние от источника до детектора.
Такой подход (урав. 6) соответствует временной области моделирования, который позволяет рассматривать влияние спеклов на сигнал с по-
мощью множителя ооб
2p(2 z - L) l
. Такой подход может быть изменен,
для того чтобы использовать его для спектрального моделирования домена, так как мы рассматриваем вес и длину оптического пути пакета, который может быть использован для регистрации сигнала [15]. В дополнение к уровню отсечки (когда вес фотонного пакет упал ниже заданного порогового значения), мы используем оптическую длину пути отсечки, когда длина оптического пути превышает значение максимальной длины когерентности источника в опорном плече интерферометра. Вместо подхода с частичным отражением, в данной модели, используется принцип 'все или ничего', поскольку он позволяет повысить точность, в то время как, общая отражательная способность объекта в обоих случаях равна. Вероятность отражения на границах кубов с разными показателями преломления, определяется уравнением Френеля:
«=1
2
sin2(ai - at) + tan2(ai - at)
2 2 sin (ai +at) tan (ai +at) у
2
где аI и аt - углы отражения и преломления, соответственно.
Моделирования структурного изображения кожи человека. Для
того, чтобы обеспечить надлежащую валидацию предложенной методики мы сравним экспериментальные изображения ОКТ с соответствующими результатами моделирования (рис. 4). В ОКТ системе используется источник излучения с центральной длиной волны 1 = 840пт и полной шириной на уровне половинной амплитуды (Е"^НМ) 50 нм [16]. Десять миллионов
<
фотонов были использованы для построения каждого А-скана. Общее количество A-сканов в моделируемой изображения 400. Использовался компьютер с четырёхядерным процессором Intel i5-4670 и 8 Гб оперативной памяти (один потока для каждого ядра).
0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.0Э 0.1
0 0.02 0.04 0.06 0.08 0.1
Рис. 4. Экспериментальное ОКТ изображение человеческой кожи после компьютерной обработки. Y-ось - глубина в см, по оси Х - ширина в см.
Оттенки серого представляют интенсивность в относительных
единицах
Оптические свойства, представленные в таблице 1, были использованы для моделирования [3, 4, 17]. Оптические свойства принимаются с соответствующими масштабными коэффициентами в соответствии с распределением интенсивности в экспериментальных изображениях. Результаты моделирования с учетом и без учета влияния спеклов, представлены на рисунке 5.
Таблица 1
Оптические свойства кожи человека
Окружающая среда M (см _1) Ma (сМ _1) & n
Роговой слой 350 0,2 0,9 1,49
Верхний эпидермис 50 0,15 0,95 1,37
Эпидермис 120 0,2 0,85 1,38
Дерма 90 0,8 0,9 1,36
Нижняя дерма 120 2 0,95 1,35
Рис. 5. ОКТ изображения человеческой кожи без влияния спеклов (а)
и со спеклами (б)
Коэффициент корреляции между А-сканами экспериментального изображения после соответствующей обработки и моделируемого изображения без учета влияния спеклов равен 0.83. Распределение интенсивности моделируемых изображений соответствует экспериментальным распределениям с высокой точностью. Роговой слой отражает подавляющее большинство фотонов так как разность показателя преломления этого слоя и воздуха очень высока (таблица 1). Кроме того, из-за высокого значения коэффициента рассеяния и фактора анизотропии, существует много фотонов, которые отражается обратно в детектор из эпидермиса. Однако, из-за применения вокселного наблюдается высокая отражающая способность на границах слоев, вычисленных в направлении нормали к границам прямых вокселей в то время как в уравнении представления биологического объекта нормали к точке рассчитывается на изогнутых поверхностях [9-10]. Исходные вокселные границы всегда прямые и их плоскости ориентированны
ортогонально ХУ7-плоскостям. Сравнение отдельных экспериментальных и моделируемых А-сканов также показывает хорошее соответствие между полученными данными (рис. 6).
Рис. 6. Сравнение одного А-скана модели и одного А-скана, взятого из экспериментального изображения. Y- ось - интенсивность в логарифмической шкале оси Х - глубина в пикселях. R = 0.85
На границе верхней роговой-воздушной границе сильное отличие отражения имеет место. Это которое может быть результатом дополнительной компенсации, выполняемой измерительной ОКТ системой. Общее время, затраченное на моделирование одного А-скана приблизительно равно одной минуте, соответственно моделирования B-скана занимает около 7-ми часов.
Заключение
Выполнение моделирования с использованием экспериментальных B-сканов в сочетании с вокселным подходом позволяет осуществлять реконструкции изображений объектов, которые представляют собой структурные изображения реальных биологических образцов.
Описанный метод моделирования структурных изображений позволяет получить А- и В-сканы, которые находятся в хорошем соответствии с экспериментальными данными. Однако, точное соответствие при моделировании может быть потеряно из-за простоты описания представленных вокселов, которые не позволяют построить произвольную форму границ объекта. Другим слабым местом описанного способа является его низкая скорость вычисления. Есть несколько других подходов к повышению скорости вычислений, включая вычисления на графических процессорах (GPU). Используемые функции рассеяния позволяют значительно увеличить скорость вычислений, однако до настоящего времени не проводились исследования, которые демонстрируют высокую точность при моделировании структурных ОКТ изображений с произвольной геометрией.
Список литературы
1. Zimnyakov D.A., Tuchin V.V. Optical tomography of tissues, Kvant. electron., 2002. 32:10. P. 849-867.
2. Proskurin S.G., Potlov A.Yu., Frolov S.V. Specific features of diffuse photon migration in highly scattering media with optical properties of biological tissues, Kvant. electron., 2015. 45:6. P. 540-546.
3. Tom Lister, Philip A. Wright, Paul H. Chappell Optical properties of human skin, J. Biomed. Opt. 17(9), 2012.
4. Ding H.F. et al. Refractive indices of human skin tissues at eight wavelengths and estimated dispersion relations between 300 and 1600 nm. Phys. Med. Biol., 2006. 51, (6 ). P. 1479 -1489.
5. Proskurin S.G. Raster scanning and averaging for reducing the influence of speckles in optical coherence tomography, Kvant. electron., 2012. 42:6 P. 495-499.
6. Wang R.K. Signal degradation by multiple scattering in optical coherence tomography of dense tissue: a Monte Carlo study towards optical clearing of biotissues. Phys. Med. Biol., 2002. 47(13) P. 2281-2299.
7. Wang L., Jacques S.L., Zheng L. MCML—Monte Carlo modeling of light transport in multi-layered tissues. Comput. Methods Programs Biomed. 1995. 47(2). P. 131-146.
8. Yao G., Wang L.V. Monte Carlo simulation of an optical coherence tomography signal in homogeneous turbid media. Phys. Med. Biol., 1999. 44(9). P. 2307-2320.
9. Kirillin M.Yu., Priezzhev A.V., Myllyla R.A. Role of multiple scattering in formation of OCT skin images, Kvant. electron., 2008. 38:6 (2008). P. 570-575.
10. Kirillin M.Yu. et al. Monte Carlo simulation of optical clearing of paper in optical coherence tomography. Quantum Electron. 2006. 36(2). P. 174180.
11. Periyasamy V., Pramanik M. Importance sampling-based Monte Carlo simulation of time-domain optical coherence tomography with embedded objects. Applied Optics, 2016. 55(11). P. 2921-2929.
12. Siavash Malektaji, Ivan T. Lima, Jr., Sherif S. Sherif Monte Carlo simulation of optical coherence tomography for turbid media with arbitrary spatial distributions. J. Biomed. Opt. 2014. 19(4).
13. Fang Q. Mesh-based Monte Carlo method using fast ray-tracing in Plücker coordinates. Biomed. Opt. Express, 2010. 1(1). P. 165-175.
14. Delia Cabrera Fernández, Harry M. Salinas, Carmen A. Puliafito Automated detection of retinal layer structures on optical coherence tomography images. Opt. Express 13, 2005.
15. Alzbeta E. Hartinger, Ahhyun S. Nam, Isabel Chico-Calero, Benjamin J. Vakoc Monte Carlo modeling of angiographic optical coherence tomography. Biomed. Opt. Express, 2014. 5.
16. Proskurin S.G., Wang R.K. Visualisation of human subcutaneous blood vessels by increasing coherence probing depth, Quantum Electron., 2004. 34:12. P. 1157-1162.
17. van Gemert M.J. C. et al. Skin optics. IEEE Trans. Biomed. Eng.. 1989. 36, (12). P. 1146 -1154.
Абдулкарим Саиф Назар, асп., saif.abdulkareem@,mail.ru, Россия, Тамбов, Тамбовский государственный технический университет,
Петров Денис Алексеевич, магистрант, den 794amail. ru, Россия, Тамбов, Тамбовский государственный технический университет,
Проскурин Сергей Геннадиевич, канд. техн. наук, доц., spros atamh.ru, Россия, Тамбов, Тамбовский государственный технический университет
MONTE CARLO SIMULA TION OF STRUCTURAL OCT IMAGES OF HUMAN SKIN IN VIVO WITH THE USE OF EXPERIMENTAL B-SCANS AND DISTRIBUTION OF OPTICAL
PROPERTIES
S.N. Abdulkareem, D.A. Petrov, S.G. Proskurin
Monte Carlo simulation of structural OCT (optical coherence tomography) images using experimental images of real bioobjects is described. The biological object is considered to be a set of 3D elements allowing simulation of the object which structure cannot be described analytically. Each voxel is characterized by its index of refraction and anisotropy parameter, scattering and absorption. Experimental acquisition of the internal structure of bio-object is used to reconstruct a simulated image instead of an analytical representation of the geometry of heterogeneities' boundaries. The efficiency of the described method is verified by comparing the simulated and experimentally acquired A- and B-scans.
Key words: optical coherence tomography, coherence probing depth, human skin, Monte Carlo simulation, voxel based model.
Abdulkareem Saif Nazar, postgraduate, saif.abdulkareema mail.ru, Russia, Tambov, Tambov State Technical University,
Petrov Denis Alekseevich, undergraduate, den 794a mail. ru, Russia, Tambov, Tambov, State Technical University,
Proskurin Sergey Gennadevich, candidate of technical sciences, docent, sprosa tamb. ru, Russia, Tambov, Tambov State Technical University
