ОПРЕДЕЛЕНИЕ МОРФОЛОГИИ ЭРИТРОЦИТОВ МЕТОДОМ РЕШЕНИЯ ОБРАТНОЙ ЗАДАЧИ СВЕТОРАССЕЯНИЯ С ИСПОЛЬЗОВАНИЕМ КЛАСТЕРИЗАЦИИ

​АННОТАЦИЯ

Разработан программный комплекс для решения обратной задачи светорассеяния с использованием кластеризации. Среди большого разнообразия был выбран иерархический метод построение бинарного дерева – метод kd-дерева. Индикатрисы, представляющие из себя точки в многомерном пространстве, разделялись на кластеры по координате, разброс амплитуд в которой оказывался наибольшим. Глубина построения дерева задается вручную оператором, и зависит от размера базы данных. Установлено, что отсутствует необходимость построения дерева до узлов единичных размеров, так как дальнейшая обработка частиц требует накопления некоторого их количества. Для использованной в работе базы данных размером 1,8 млн. частиц оптимальная глубина построения составляет 18. В процессе работы для большего ускорения реализован промежуточный этап фильтрации узлов. Фильтрация производилась по параметру связи теоретического и экспериментального сигнала, названного в работе коэффициент α.

Ключевые слова: кластеризация, индикатриса светорассеяния, сканирующий проточный цитометр, обратная задача светорассеяния.

Исследование морфологии здоровых и патологически измененных эритроцитов может помочь осуществлять более раннюю и точную диагностику болезней, а также исследовать взаимосвязь между известными болезнями и параметрами эритроцитов. В научной среде проведение экспериментов, как правило, не так сильно ограничено по длительности и более требовательно к точности осуществляемых исследований. Из этого следует, что не все методы, используемые в науке, пригодны для диагностики. В настоящее время в клиниках чаще всего используют микроскопию [1, с. 7; 2, с. 405; 3, с. 623; 4, с. 1040], кондуктометрический метод и проточную цитометрию.

Сканирующая проточная цитометрия – метод исследования частиц путем анализа рассеянного ими света. Проба, содержащая исследуемые частицы, прогоняется через капилляр, освещаемый лазером, частицы рассеивают свет лазера, рассеянный свет при помощи системы зеркал доставляется в фотоэлектронные умножители. Далее измеренный сигнал обрабатывается для получения индикатрисы светорассеяния – зависимости интенсивности светорассеяния от угла рассеяния. Индикатрисы содержат в себе всю необходимую информацию о светорассеивающих свойствах частицы.

Актуальный метод исследования при помощи сканирующей проточной цитометрии на базе лаборатории цитометрии и биокинетики (ЦиБ) СО РАН позволяет измерять вплоть до 48 параметров эритроцитов, что значительно больше, чем результаты измерений на стандартном цитометре. Однако для практических целей существенным недостатком данной технологии является большое время обработки данных. Обработка результатов измерения пробы для одного пациента в настоящее время занимает более часа, что неприемлемо в клинических исследованиях. По этой причине было решено исследовать возможности ускорения алгоритма для дальнейшего его внедрения в клиническую систему.

Для определения параметров измеренных частиц необходимо провести решение обратной задачи светорассеяния. Экспериментально измеренные на цитометре индикатрисы сравниваются с предварительно насчитанной теоретической базой данных частиц с известными параметрами, и в результате сравнения параметры наиболее близкой (с точки зрения нормы) частицы из базы данных принимаются с некоторой вычисляемой погрешностью за параметры исследуемой частицы. Таким образом, проводя сравнения, находятся параметры для всех измеренных частицы в пробе.

В работе реализован метод ускорения решения обратной задачи светорассеяния, основанный на кластеризации используемой базы данных методом kd-дерева, а также введена дополнительная фильтрация узлов построенного дерева.

Метод kd-дерева – алгоритм кластеризации, осуществляющий разбиение данных на двоичное дерево, в котором каждый узел является k-мерной точкой. Каждый нелистовой узел можно рассматривать как гиперплоскость, которая делит пространство на две части, известные как полупространства. Точки слева от этой гиперплоскости представлены левым поддеревом этого узла, а точки справа от гиперплоскости - правым поддеревом.

Пусть имеется некоторый набор данных, состоящий из p точек N‑мерного пространства V₀ = {v_1, …, v_p}, v_i = (xⁱ₁, …, xⁱ_N)  R^N. Для того, чтобы разделить текущее дерево на 2 листа, найдем такую координату k в пространстве R^N, для которой:

                                                                              (1)

Далее, отсортируем все точки в V по возрастанию значения в координате k, найдем медианное значение v_med и сформируем два новых кластера: в левый кластер V₁ попадут все точки отсортированного кластера от начала до v_med не включительно, в правый кластер V₂ попадут все точки от v_med до последней точки.

Используя вышеописанный алгоритм, можно построить дерево нужной глубины. В стандартных случаях, дерево подразумевается простроенным “до конца”, т.е. кластеры самого нижнего уровня будут содержать в себе по одному элементу исходного набора данных. Однако в случае решения обратной задачи светорассеяния такой подход хоть и возможен к реализации, но оказался менее эффективным по сравнению с реализованным в данной работе.

Далее построенное дерево будет использовано для решения вариации задачи поиска ближайшего соседа – поиск ближайшего диапазона элементов. Важнейшим параметром для решения такой задачи является мера расстояния.

Амплитуда экспериментального сигнала и смоделированного сигнала отличается на порядки. Следовательно, для корректного их сравнения необходимо ввести некий коэффициент, названный в данной работе коэффициентом α. По своей физической сути, коэффициент α – это коэффициент экспериментальной установки, зависящий от многих факторов.

Математический вывод этого коэффициента основан на решении уравнения:

,                                                             (2)

где exp_i и th_i – экспериментальная и теоретическая индикатрисы, соответственно, ω_i – инструментальный весовой коэффициент.

 Для решения этого уравнения необходимо взять производную по α приравнять к нулю:

                                                  (3)

Отсюда коэффициент α равен:

                                                                                    (4)

Так как этот коэффициент связан с настройками установки, перед проведением эксперимента при помощи модельных частиц (полистирольных микросфер с хорошо известными параметрами), рассчитывается диапазон возможных коэффициентов α для данной экспериментальной пробы. Знание об этом диапазоне равнозначно знанию того, какие амплитуды экспериментальных индикатрис будут приемлемы с физической точки зрения.

Алгоритм реализован при помощи среды разработки Labview. Работа реализованной программы и все дальнейшие результаты были получены на стационарном ПК со следующими характеристиками: процессор Xeon Е5-2689 тактовой частотой 2,6 Ггц, 256 гб оперативной памяти.

Все вычисления проводились с использованием базы данных эритроцитов размером 1,8 млн. индикатрис, набор экспериментальных индикатрис составлял 2000 частиц. В настоящее время цитометр устроен таким образом, что за один экспериментальный цикл измеряется порядка 2000 частиц, так как такое количество можно считать статистически достоверным.

Первые опыты с кластеризацией показали, что выигрыш по времени составляет порядка 10-15% по сравнению с актуальным методом лаборатории. Положение уровня дополнительной фильтрации непосредственно влияет на точность и скорость работы всего алгоритма. Более высокое положение узла в дереве свидетельствует о том, что индикатрисы данного узла претерпели меньше процедур разделения, чем все нижележащие узлы. Следовательно, формы индикатрис в более высоких узлах внутри самих этих узлов отличаются сильнее, чем в нижележащих. А значит, средняя индикатриса в высоких узлах, сильнее отличается от представителей такого узла. Из этого следует, что чем раньше применяется фильтрация, тем больше становится вероятность того, что потенциально “хорошие” индикатрисы могут быть отсеяны, так как попали в кластер с “плохой” средней индикатрисой. Таким образом, необходимо решить задачу оптимизации: определить такую глубину фильтрации, чтобы отношение точности и скорости работы было удовлетворительным.

В качестве эталона работы реализованного алгоритма был принят актуальный лабораторный метод, как более точный и проверенный. Таким образом, ожидалось, что при увеличении глубины фильтрации результаты работы алгоритма с применением кластеризации будут стремиться к результатам лабораторного алгоритма.

Для тщательного сравнения результаты были представлены в следующем виде:

,                                                                         (5)

где ,  – экспериментальный и теоретический векторы параметров для i-ой индикатрисы, соответственно. Выбранное представление показывает отклонение средних параметров, вычисленных методом кластеризации, от эталонного значения. Результаты представлены в таблице 2 и на рисунке 8. Как видно, наибольшее отклонение от эталона наблюдается при измерении параметра толщины перетяжки, и составляет порядка 30%. Такое высокое отклонение связано с тем, что эталонный метод изначально плохо определяет данный параметр, посему предлагается исключить его из рассмотрения.

Таблица 1.

Показатели отклонений (%) основных параметров модели при различных показателях глубины фильтрации: d – диаметр, w – толщина перетяжки, h – наибольшая толщина, n – показатель преломления, V – объем эритроцита.

Рисунок 1. Зависимость отклонения параметров эритроцитов от глубины фильтрации: d – диаметр, h – наибольшая толщина, n – показатель преломления, V – объем

Из сравнения отклонений видно, что при увеличении глубины фильтрации значения параметров при анализе кластерным методом сходятся к эталонным. Таким образом, для определения оптимального или субоптимального уровня глубины, необходимо выбрать, какое время работы алгоритма будет удовлетворять потребностям в клинической практике. На рисунке 9 и в таблице 3 представлена зависимость среднего времени обработки одной индикатрисы от глубины фильтрации.

Рисунок 2. Зависимость среднего времени обработки индикатрисы от глубины фильтрации по параметру альфа

Из результатов видно, что при увеличении глубины фильтрации, начиная с уровня 14, происходит резкое повышение времени работы алгоритма. Отсюда следует, что субоптимальным решением будет уровень фильтрации 14, поскольку увеличение времени работы, по сравнению с предыдущими уровнями, незначительное, однако такое решение лежит ближе к эталону. Субоптимаьным оно является из-за того, что решение все же отличается от эталона, а не совпадает с ним, однако выигрыш по времени в случае такого анализа, преобладает над точностью решения.

Таблица 2.

Время обработки одной индикатрисы при различных уровнях фильтрации

Список литературы:

1. Bransky A. et al. Correlation between erythrocytes deformability and size: a study using a microchannel based cell analyzer // Microvasc. Res. 2007. Vol. 73, № 1. P. 7–13.
2. Canham P.B., Burton A.C. Distribution of Size and Shape in Populations of Normal Human Red Cells // Circ. Res. 1968. Vol. 22, № 3. P. 405–422.
3. Gifford S.C. et al. Parallel Microchannel-Based Measurements of Individual Erythrocyte Areas and Volumes // Biophys. J. 2003. Vol. 84, № 1. P. 623–633.
4. Tomaiuolo G. et al. Comparison of two flow-based imaging methods to measure individual red blood cell area and volume // Cytom. Part J. Int. Soc. Anal. Cytol. 2012. Vol. 81, № 12. P. 1040–1047.