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

В последнее время информационные и цифровые технологии находят все большее применение в разнообразных видах деятельности человека. Одним из эффективных применений является анализ свойств оптически проницаемых сред. В задачах контроля качества и управления технологическим процессом гомогенизации молочных продуктов может использоваться микроскопия проб молочных продуктов для анализа степени гомогенизации, распределения жировых шариков по массе [1]. Для автоматизации обработки информации используются микроскопы, соответствующие цифровые камеры, обеспечивающие регистрацию изображений в заданные промежутки времени, устройства захвата изображений или видеопотока и управления цифровой камерой, а также программные средства для интерпретации полученной информации в виде, удобном для лица, принимающего решение (ЛПР).

Несмотря на наличие средств получения цифровых микроснимков проб молочных продуктов, остается актуальной задача разработки специального программного обеспечения (СПО), реализующего весь цикл обработки информации применительно к молочным продуктам. Для разработки СПО можно использовать имеющиеся библиотеки и среды обработки цифровых изображений, которые необходимо скомпоновать в соответствии с заданным алгоритмом [2]. Одной из распространенных и удобных сред, позволяющих создавать гибкие алгоритмы обработки графической информации является среда ImageJ, в которой реализовано большое количество типовых алгоритмов обработки и анализа цифровых изображений [3]. Кроме того, использование механизма макросов на языке Java позволяет создавать специфические алгоритмы обработки цифровых изображений в соответствии с заданной задачей, импортировать данные разных форматов, а также экспортировать данные в продукты семейства Microsoft Office и Statistica, делает этот продукт удобным для отработки алгоритмов анализа цифровых изображений в задачах анализа качества молочных продуктов.

Рассмотрим использование среды ImageJ для задачи исследования изменения распределения жировых шариков молока в процессах гомогенизации молока.

Процесс гомогенизации молока и жидких молочных продуктов служит для повышения дисперсности в них жировой фазы, что позволяет исключить отстаивание жира во время хранения молока, развитие окислительных процессов, дестабилизацию и подсбивание при интенсив- ном перемешивании и транспортировании. Под степенью гомогенизации понимают средний размер жировых шариков, степень устойчивой во времени однородности (гомогенности) в смысле отсутствия микроконцентрационных неоднородностей, образующихся при смешивании взаимно-нерастворимых веществ [1].

Изменение распределения жировых шариков осуществлялось гомогенизацией исходного образца молочного продукта (фермерского негомогенизированного молока) на лабораторной установке гомогенизации научнопроизводственной лаборатории группы компаний «ЭФКО», позволяющей осуществлять процесс гомогенизации молочных продуктов с разным давлением гомогенизации. Для исследований использовались давления 30, 60, 90, 120, 150, 180 МПа. Полученные после гомогенизации пробы исследовались с помощью цифрового микроскопа с увеличениями 100, 400, 600 крат. Параллельно осуществлялось по 2 снимка, ввиду чего общее количество снимков составило 42 штук (таблица 1). Пример микроснимка пробы молока с увеличением в 100 раз представлен на рисунке 1.

Т а б л и ц а 1

Основные характеристики микроснимков

№ п/п	Давление гомогенизации, МПа	Увеличение, крат
1	Без гомогенизации	100
2	Без гомогенизации	400
3	Без гомогенизации	600
4	30	100
5	30	400
6	30	600
7	60	100
8	60	400
9	60	600
10	90	100
11	90	400
12	90	600
13	120	100
14	120	400
15	120	600
16	150	100
17	150	400
18	150	600
19	180	100
20	180	400
21	180	600

Для калибровки изображений и перехода к физическим единицам измерения осуществлялась микросъемка эталонного объекта (рисунок 1).

Поскольку цвет пробы не представляет интереса для дальнейших исследований, для упрощения дальнейших преобразований необ- ходимо осуществить преобразование в оттенки серого с желаемой шкалой, например Yy = 0.3Ry + 0.59Gy + 0.11 B . Осуществим преобразование в режим 32 bit следующей командой: run("32-bit").

Дальнейшим этапом является приведение цифрового изображения к физическим единицам измерения с использованием шкалирования на основе цифрового изображения измерительной линейки (рисунок 2) с одинаковым коэффициентом увеличения. Для этого используется инструмент MakeLine с параметрами, которые автоматически задаются при создании линии с длиной, равной одному делению эталонной шкалы. Соответствующий скрипт макроса на языке Java имеет вид:

makeLine(x1, y1, x2, y2);

run("Set Scale...", "distance=26.019 known=0.000001 pixel=1 unit=m");

где x1, y1, x2, y2 – координаты точек начала и конца отрезка, соответствующего одному делению шкалы (рисунок 2), параметр distance – длина отрезка в пикселях, known – длина отрезка в единицах unit.

Например, одному делению шкалы на рисунке 1 при увеличении 100 раз соответствует 26 пикселей. В дальнейшем все снимки проходят процедуру установки масштаба Ana-lyze->SetScale… с параметрами Distanceinpixels = 26, Pixelaspectratio = 0.000001 и Unitoflength = m.

Далее необходимо избавиться от помех на фоне, так как некоторые мелкие объекты в дальнейшем могут ошибочно распознаться как жировые шарики. Для этого можно использовать функцию Process->SubtractBackground… с параметром Rollingballradius =170.

run("SubtractBackground...", "rolling=150");

Рисунок 2. Снимок калибровочной линейки со 100 кратным увеличением

Альтернативным вариантом является использование полосового фильтра (bandpassfilter) на основе быстрого преобразования Фурье, чтобы отфильтровать объекты в заданном диапазоне их размеров (например, 5-600px). Алгоритм осуществляет переход в частотную область с использованием быстрого преобразования Фурье, осуществляет вырезку необходимого частотного диапазона, а затем выполняется обратное преобразование[4].

Изображение f ( x, y ) измеряется с шагом A x , A y , измеренные значения f_jx^ содержат полезный сигнал и помеху [2]

где 1 < jx < Nx ,1 < jy < Ny , Nx. Ny - размер изображения в единицах элементов размером Ax xAy (может быть один пиксель). Далее осуществляется вычисление коэффициентов дискретного преобразования Фурье (ДПФ):

\       1     N x N y -        - 2 7V ⁽ j-l- ⁺ j y l y )

F ( l x . ' . ) = — Z Z f,,,-e^NN

N_XN₂ j_x = o j_y = o

Каждому индексу l , l соответствует ча-

Рисунок 1. Снимок пробы молока со 100-кратным увеличением

^стота to x =^A _M, l x . to y =^A _% l y .

Далее из полученного спектра изображения удаляются определённые участки за пределами полосы пропускания ^п <  to <  to , где ^ , to - нижний и верхний пределы полосы пропускания, to — частоты, попадающие в диапазон полосы пропускания. Для этого осуществляется операция умножения спектра изображения на функцию, определяющую характер «вырезания» определённых частот

где W ( l_x. l_y ) функция, реализующая удаление части спектра. В простейшем случае закон «вырезания» выглядит следующим образом:

0, V w >  Ю нп

L (to) = П, V^ > to > Юнп

0, ^{V w} <  « вп

Этому преобразованию соответствует оператор:

run("BandpassFilter...", "filter_large=600 filter_small=5 suppress=Nonetolerance=5 au-toscalesaturate").

Результирующее изображение представлено на рисунке 3.

Следующим этапом идет процедура получения двухбитного изображения на основе пороговой классификации пикселей в соответ- ствии с правилом:

T = { P ^e Df ( P ) < h } ,         ⁽⁵⁾

где p – характеристика «высоты» пикселя; h – величина порога.

Для этого используются процедуры:

setThreshold(1, 159); run("ConverttoMask").

Следует отметить, что выбор пороговых значений в сильной степени определяет качество анализа частиц в дальнейших процедурах, т.к. в результате этой операции происходит как «склеивание» рядом лежащих частиц в одну, так и «незамыкание» колец для светлых частиц, что приводит к погрешностям в дальнейшем анализе. Такие проблемы вызваны как неравномерностью засветки по всей площади снимка, так и ввиду того, что проба молока имеет свою толщину и жировые шарики находятся на разном расстоянии от точки фокуса и поэтому имеют разную резкость и яр-кость.Результирующее изображение представлено на рисунке 4. Для автоматического определения порогов можно использовать процедуру setAutoThreshold("Default").

Рисунок 3. Вид микроснимка с настроенным диапазоном полосы пропускания

Рисунок 4. Микроснимок в двухбитном виде

Поскольку жировые шарики частично пропускают свет, после обработки получаются кольцевые структуры. Для избавления от кольцевых структур можно использовать операцию заливки колец (рисунок 5) run("FillHoles").

Наличие «слипшихся» элементов является причиной близкого расположения жировых шариков друг к другу и их частичного перекрытия. В дальнейшем изображение после операции FillHoles малопригодно для последующей обработки, т.к. «слипшиеся элементы» в ряде случаев будут интерпретироваться как одно целое и необходима процедура, снижающая погрешности расчетов, вызванные «слипанием». Для этих целей можно использовать алгоритм преобразования графических изображений Watershed, осуществляющий декомпозицию отдельных элементов изображения по принципам геометрических характеристик их пикселей, аналогичных высоте в рельефе земли и моделирующий «водораздел» [5]. Для бинарных изображений алгоритм декомпозиции основан на вычислении расстояний каждой точки изображения из множества черных пикселей b e B до ближайшей точки из множества белых пикселей w e W . Расстояние в простейшем случае может быть Евклидовым: D eu ( [ i i , j i ] , [ i 2 , j 2 ] ) = V( i i - i 2 ) 2 + ( j i — j 2 ) 2 . (6)

Преобразование осуществляется процедурой run("Watershed"). Результат преобразования представлен на рисунке 6.

Рисунок 5. Изображение с залитыми кольцами

Рисунок 7. Изображение в виде контуров

Рисунок 6. Изображение, преобразованное по принципу геометрических характеристик

Полученная таблица содержит информацию о площади обнаруженного объекта в выбранной на этапе шкалирования размерности (пиксель2, или мкм2), среднем, минимальном и максимальном значении выборки, а также порядковый номер объекта. Формат Excel сохранения таблицы совместим с пакетом Statistica,

что позволяет использовать всю вычислитель-

ную мощь этого пакета для дальнейшего ана-

лиза полученных данных.

Для дальнейшего анализа также целесообразно осуществить переход от единиц пло-

щади к радиусу объекта

r o6 = ^

и, зная

плотность жира, к приблизительному значе-

нию массы m

об    ^{п Пр} ж ^Гоб

Проведенные операции подготавливают изображение к проведению основной операции, которая осуществляет поиск и анализ частиц на изображении в соответствии с заданными параметрами AnalyzeParticles. При использовании этой функции можно задать диапазон размеров площадей частиц, которые будут анализироваться в пикселях в квадрате (Size, pixel^2), причем максимальный диапазон может быть от 0 до значения Infinity. Параметр Circularity отвечает за анализ близости найденных частиц к форме круга, значение 1 соответствует кругу, а ноль - полосе единичной толщины. Промежуточные значения позволяют настроить алгоритм поиска на обнаружение фигур, близких по форме к кругу (например, Circularity=(0.6÷1)).

Пример выходного графического изображения с опцией Outlines (контуры) представлен на рисунке 7. Также в результате работы процедуры AnalyzeParticles формируется таблица,содержащая перечень свойств всех обнаруженных объектов.

Эти расчеты можно выполнить как в пакете Statistica, так и в системе Excel.

Обработка преобразованных данных в пакете Statistica имеет множество возможно-

стей. Основными представляющими интерес с точки зрения свойств молочных продуктов и влияния технологического режима их обработки (в данном случае процесса гомогенизации) являются простейшие описательные статистики (descriptivestatistics). Выделением соответствующего столбца данных и выбора меню Statistics рассчитываются такие статистические характеристики исследуемой выборки

как среднее значение, размах, дисперсия, коэффициенты асимметрии и эксцесса.

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

показывающие отклонение от нормального распределения и «ящичковые»-диаграммы.

Рисунок 8. Гистограмма распределения жировых шариков в анализируемом образце по радиусу

В заключение можно отметить, что в сравнении с результатами визуального анализа, при обработке «потеряно» порядка 7% шариков, что можно считать приемлемым результатом, основной проблемой подхода является наличие незамкнутых кольцевых областей, которые не отрабатываются процедурой fillholes. Улучшить результаты можно многократным применением и комбинированием фильтров.