Сенсорные системы, 2021, T. 35, № 4, стр. 340-347

Генерирование синтетических пористых изображений для аугментации данных с целью тренировки алгоритмов машинного обучения

А. В. Хафизов 1*, М. В. Григорьев 2

1 ФНИЦ Кристаллография и фотоника РАН
119333 Москва, Ленинский проспект, 59, Россия

2 ИПТМ РАН
142432 Черноголовка, ул. Академика Осипьяна, 6, Россия

* E-mail: ankhafizov1998@yandex.ru

Поступила в редакцию 09.06.2021
После доработки 16.07.2021
Принята к публикации 03.08.2021

Полный текст (PDF)

Аннотация

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

Ключевые слова: компьютерная томография, трансферное обучение, аугментация, синтетические данные

ВВЕДЕНИЕ

Последнее время проводится множество исследований, посвященных изучению пористых структур методом компьютерной абсорбционной томографии (КТ). Эти исследования являются перспективными, потому что пористые материалы имеют широкое применение в науке и индустрии. Например, пористые металлокерамические мембраны используются в очистке жидкостей (Uvarov et al., 2014), полимерные материалы – для изготовления биоразлагаемых имплантов (Borovikov et al., 2018).

Томография может использоваться для контроля геометрических параметров пор, которые отвечают за функциональные свойства рассмотренных материалов. Этот метод исследования является не деструктивным, т.е. в отличие от множества других более консервативных способов исследуемый образец не подвергается полному или частичному разрушению. Более того, при помощи КТ можно получить цифровое 3D-изображение образца и использовать его для моделирования физических явлений, например, протекания жидкости сквозь него (Shah et al., 2016).

Компьютерная томография имеет недостаток – значительное количество времени, необходимое для получения данных, особенно при экспериментах на лабораторных источниках излучения. Так, для съемки одного образца может понадобиться время, порядка нескольких часов. Сократить время можно, если уменьшить количество радоновских проекций или время экспозиции для каждой проекции (Buzug, 2008). Однако сокращение времени таким образом приводит к увеличению уровня шума на реконструированном КТ изображении, что может негативно повлиять на его анализ, в частности, на сегментацию и бинаризацию – определение принадлежности каждого вокселя к классу пор или материала. Тем не менее можно попробовать уменьшить шум при предобработке экспериментальных данных, и для этого существует много техник фильтрации (Usanov et al., 2017; Chukalina et al., 2021). Приведенные методы принадлежат к группе так называемых консервативных методов компьютерного зрения, которые были известны еще в 1990–2010 годах. Они являются относительно простыми, но зачастую требуют значительных вычислительных ресурсов и/или являются параметрическими, а значит, требуют от исследователя опыта применения и умения разбираться в особенностях этих методов.

На данный момент существуют более современные способы обработки и анализа изображений, основанных на машинном обучении (ML), которые показывают лучшие результаты с точки зрения качества и затрат времени на бинаризацию. Часть из них основана на классическом ML, например, на применении логистической регрессии для бинаризации пористых структур (Хафизов и др., 2020), другая часть является более сложными c использованием экстракции признаками вручную (Jaccard et al., 2017) или сверточных нейронных сетей (CNN) (Ronneberger et al., 2015; Badrinarayanan et al., 2017). Стоит отметить, что CNN, хотя и являются быстродействующими и непараметрическими, но все же обладают низкой устойчивостью и не всегда работают “из коробки” (Brudfors et al., 2021). Для лучшего качества бинаризации можно переобучить CNN на данных, которые требуется бинаризовать. Аналогично можно дообучать более простые ML-алгоритмы (Хафизов и др., 2020; Jaccard et al., 2017). Таким образом, все предложенные методы требуют предварительного обучения перед их использованием.

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

Тренировка CNN может быть реализована с использованием идеи трансферного обучения, когда при обучении изменяются веса только последних слоев, а не всей ML-модели. Этот метод является перспективным, так как он позволяет наименьшими затратами времени, и используя небольшой датасет для обучения, получить рабочую нейросеть для частного случая данных. Тем не менее относительно много пар изображений все еще необходимо для обучения, что тяжело получить вручную.

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

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

ОПИСАНИЕ ГЕНЕРАТОРА ФАНТОМА

В качестве генератора пористых структур мы использовали модифицированный алгоритм blobs из пакета PoreSpy (Gostick et al., 2019). Отличием нашей версии является то, что вместо безразмерного нормированного параметра blobiness, определяющего характерные размеры пор/материала, мы ввели параметр, напрямую определяющий такие размеры в пикселях. На вход генератор принимает три параметра: размеры генерируемого изображения, желаемую пористость p (отношение объема пор к общему размеру образца) и σ (среднеквадратичное отклонение сглаживающей гауссианы), которая определяет характерные размеры. Схема алгоритма генерации такова:

1. Создание массива заданного размера и заполнение его случайными значениями из интервала [0.0, 1.0) с равномерным вероятностным распределением.

2. Свертка массива из п. 1 с гауссианой со средним квадратичным отклонением (СКО), равным σ.

3. Гистограммная эквализация (Gonzalez, Woods, 2002) массива, полученного в п. 2.

4. Пороговая бинаризация с порогом, равным пористости p.

Пример фантомов с разными значениями параметров σ и p показан на рис. 1.

Рис. 1.

Примеры фантомов размера 1000 × 1000 вокселей с параметрами.

а – {p = 0.5, σ = 30}; б – {p = 0.5, σ = 10}; в – {p = 0.2, σ = 10}. Поры – черные.

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

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

Для исследования устойчивости алгоритмов к шуму эксперименты проводили на зашумленных фантомах. Для зашумления изображений использовали метод “соль и перец” (Gonzalez, Woods, 2002) как наиболее простой и схожий с шумом на экспериментальном изображении после бинаризации. В качестве верхнего порогового параметра зашумления, определяющего количество шума, выбраны значения 0 (т.е. шум отсутствует), 0.01, 0.05, 0.1, 0.2, 0.3. За нижний пороговый параметр была выбрана величина 1 – u. На рис. 2 представлены фрагменты 2D-фантомов для демонстрации уровня шума.

Рис. 2.

Фантомы, зашумленные методом “соль и перец”.

аu = 0; бu = 0.01; вu = 0.05; гu = 0.1; дu = 0.2; еu = 0.3.

Исходя из нашего опыта, если перед бинаризацией экспериментального изображения была осуществлена его предварительная фильтрация, то шум редко превышает уровень, представленный на рис. 2, г, т.е. u = 0.1.

МЕТОДИКА ОПРЕДЕЛЕНИЯ ПОРИСТОСТИ

Пористость (отношение объема пустого пространства к объему образца) в данной работе определяется по одной из методик, описанных нами в статье (Григорьев и др., 2021). В этой публикации были проанализированы несколько вариантов определения пористости, которые были ранжированы по точности работы и времени выполнения. В настоящей работе мы используем способ, который имеет оптимальное соотношение точность–время. Его сущность заключается в последовательном применении операций локальной усредняющей фильтрации, пороговой бинаризации Оцу (Kurita et al., 1992) и бинарных морфологических операций закрытия, а затем открытия ядром тетрагональной структуры. В результате применения описанного алгоритма пористость полученного бинарного изображения достаточно хорошо совпадает с реальной пористостью, несмотря на то, что могут быть размыты границы пор, присутствовать артефакты в виде “висячих камней” и др. (Chukalina et al., 2021).

МЕТОДИКА ОПРЕДЕЛЕНИЯ СКО ГАУССОВА ФИЛЬТРА

Методика определения СКО основывается на вейвлет-анализе изображений.

Вейвлет преобразованием (ВП) одномерного пространственного сигнала S(t), где t – пространственная координата, называется его представление в виде обобщенного ряда или интеграла Фурье по системе базисных функций вида (Gonzalez, Woods, 2002):

${{\psi }_{{a,b}}} = \frac{1}{{\sqrt a }}\psi \left( {\frac{{t - b}}{a}} \right),$
которые сконструированы из исходного вейвлета Ѱ(t) за счет сдвига по времени на b и изменения масштаба в a раз.

Непрерывное прямое вейвлет преобразование функции S(t) можно представить следующим образом:

${{W}_{s}}(a,b) = \frac{1}{{\sqrt a }}\int\limits_{ - \infty }^{ + \infty } {S(t)} \left( {\frac{{t - b}}{a}} \right)d$
и обратное

$S(t) = \frac{1}{{{{C}_{\psi }}}}\int\limits_{ - \infty }^{ + \infty } {\int\limits_{ - \infty }^{ + \infty } {{{W}_{s}}(a,b)} {{\psi }_{{a,b}}}(t)} \frac{{dadb}}{{{{a}^{2}}}}.$

Величина Ws(a, b) является вейвлет-спектром сигнала S(t) и зависит от масштабирующего параметра a и сдвигового b. Спектр с маленькой величиной a отражает более мелкие детали сигнала S(t), тогда как большая величина a отвечает за основные детали.

Так как непрерывное прямое ВП представляет собой разложение сигнала по вейвлет-функциям, то для наилучшей аппроксимации в разложении использовался вейвлет Рикера (Shirani et al., 2021):

$\psi (t) = \frac{1}{{\sqrt a }}\exp \left( { - \frac{{0.5{{t}^{2}}}}{{{{a}^{2}}}}} \right)({{t}^{2}} - 1).$

Вид функции выбранного вейвлета представлен на рис. 3.

Рис. 3.

Вейвлет Риккера (Shirani, 2021).

Любое КТ изображение I пористого образца можно представить как n рядов вокселей $I = \{ {{S}_{k}}(t)\} _{{k = 1}}^{n}$. Как правило, каждое отдельное Sk(t) имеет вид, представленный на рис. 4, a. Каждая яма соответствует малому значению коэффициента рентгеновского поглощения, т.е. поре, а плато – материалу.

Рис. 4.

Непрерывное вейвлет преобразование.

а – исходный зашумленный профиль КТ-изображения, б – визуализация матрицы коэффициентов после ВП сигнала (а) с вейвлетом Риккера.

Разложение сигнала из рис. 4, а при помощи вейвлета Рикера показано на рис. 4, б. Если масштабный параметр a близок к 0, то Ws(a, b) несет в себе информацию о малых возмущениях, т.е. о шуме. При больших a спектр отражает крупные детали сигнала – поры и материал.

Идея нашего подхода состоит в том, что для любого сигнала S(t) можно найти такое a0, что только одна спектральная линия Ws(a0, b) будет наилучшим образом соответствовать структуре крупных деталей этого сигнала. Тогда соответствующая a0 будет коррелировать со средним размером пор на сигнале.

Для определения наилучшего Ws(a0, b) мы пользовались коэффициентом Пирсона (Benesty, 2009):

$r(S,{{W}_{s}}({{a}_{0}},b)) = \frac{{\operatorname{cov} (S,{{W}_{s}}({{a}_{0}},b))}}{{{{\Theta }_{S}}{{\Theta }_{{{{W}_{s}}}}}({{a}_{0}},b)}},$
где Θ(.) – среднеквадратичное отклонение величины, a cov(S, Ws(a0, b)) – ковариация между S(t) и Ws(a0, b). Величина r, равная 1, означает полное совпадение S и Ws(a0, b), поэтому оптимальной величиной параметра a0 спектра сигнала будем считать такую, что:

$r(S,W({{a}_{0}},b)) \to {{\max }_{{(6)}}}.$

Из анализа зависимости a0 и σ при разных пористостях было эмпирически выявлено, что

$\left\{ \begin{gathered} d = \min (p,1 - p) \hfill \\ \sigma = ( - 3.29{{d}^{3}} + 4.31{{d}^{2}} - 1.9d + 0.9){{a}_{0}}. \hfill \\ \end{gathered} \right.$

Таким образом, предлагаемый алгоритм по определению σ(I) вейвлет-анализом на изображении I состоит из следующих шагов:

– Раскладываем I на одномерные сигналы {Sk(t)}k.

– Для каждого сигнала ищем a0, пользуясь критерием (6).

– Для каждого сигнала ищем σk и получаем множество {σ k}.

– Ищем медиану по {σ k}. Получаем σ(I).

Для 2D- и 3D-образцов значения σ можно найти вдоль каждой из декартовых осей изображения образца. При отсутствии анизотропии формы и размеров пор все значения должны быть примерно равны друг другу.

ТЕСТИРОВАНИЕ РАБОТЫ АЛГОРИТМА НА ФАНТОМАХ

Определение качества работы вейвлет поиска σ производилось на одномерных фантомах, о которых уже было упомянуто выше. Были сгенерированы фантомы с различными значениями триплетов параметров {p, σ, u}, и для каждого из них был использован предложенный алгоритм для поиска σ. Результаты численных экспериментов представлены на рис. 5: приведены графики для пористостей 0.1, 0.3, 0.5. Следует отметить, что для пористостей, больших 0.5, бинарное изображение можно представить как инвертированную картинку меньших пористостей, таким образом, например, для пористостей 0.9 и 0.7 результат будет аналогичен таковому, полученному при пористостях 0.1 и 0.3 соответственно.

Рис. 5.

Результаты работы алгоритма для определения СКО σ генератора на фантомах разных значений пористостей p и шумов u.

а – пористость 10%; б – пористость 30%; в – пористость 50%.

Как видно из сравнений рис. 5, а–в, алгоритм наиболее устойчив к шуму при пористости 0.5. Напротив, при пористости 0.1 возникают отклонения предсказанных значений от заданных, которые особенно заметны при высоких значениях шума. Это связано с тем, что при пористостях, близких к 0 (или 1), поры (или материал) мало отличаются от шума по размерам, что и вносит искажения в итоговый результат.

Для количественной оценки качества работы алгоритма использовалась метрика MAPE (Myttenaere et al., 2016), которая представляет собой среднее относительных отклонений теоретических значений компонентов n-мерной векторной величины y от экспериментально полученных значений x:

${\text{MAPE}}(x,y) = \frac{1}{n}\sum\limits_{k = 1}^n {\frac{{{\text{|}}{{x}_{k}} - {{y}_{k}}{\text{|}}}}{{{{y}_{k}}}}} .$

На рис. 6 показана зависимость значений MAPE (σзаданная, σрассчитанная) для заданного и рассчитанного параметра СКО фантома. Как видно, при шуме 0.1 и ниже наблюдается достаточно высокая точность работы алгоритма.

Рис. 6.

Количественная оценка определения СКО фантома методом вейвлет-анализа.

РЕЗУЛЬТАТЫ ПРИМЕНЕНИЯ АЛГОРИТМА НА ЭКСПЕРИМЕНТАЛЬНЫХ ИЗОБРАЖЕНИЯХ

Методика определения параметров для генерации пористых фантомов была также испытана на экспериментальных КТ-данных: реконструированных изображениях металлокерамической мембраны (на SiC основе) и полилактидных полимеров PDL разных пористостей (рис. 7, a–c). Как видно, все три образца имеют разные значения пористости и размеров пор.

Рис. 7.

Экспериментально полученные КТ-изображения (а–в) и соответствующие им синтезированные фантомы (г–е).

Пористость образцов составила приблизительно 0.48 для PDL-05(1), 0.34 для SiC мембраны и 0.52 для PDL-05(2). Параметры σ были рассчитаны по предложенному алгоритму вейвлет-анализа изображений и составили 4.51, 5.94 и 22.07 для PDL-05(1), SiC мембраны и PDL-05(2) соответственно.

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

ВЫВОДЫ

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

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

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

Поиск σ осуществлялся при помощи вейвлет преобразования. Алгоритм также показал свою высокую эффективность на численных экспериментах с одномерными протяженными фантомами, где уменьшено влияние граничных условий. Было показано, что алгоритм приемлемо устойчив к шуму на изображении.

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

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

  1. Григорьев М.В., Назиров И.В., Могилевский Е.И., Хафизов А.В., Чукалина М.В. О проблемах при работе с микротомографическими изображениями пористых структур, используемыми для моделирования процессов протекания. Труды института системного анализа российской академии наук. 2021. Т. 71. № 1. https://doi.org/10.14357/20790279210110

  2. Хафизов А.В., Григорьев М.В. Метод логистической регрессии для бинаризации томографических изображений пористых объектов. Современные методы анализа дифракционных данных и актуальные проблемы рентгеновской оптики: 7-я международная школа-семинар молодых ученых. 2020. С. 162.

  3. Badrinarayanan V., Kendall A., Cipolla R. Segnet: a deep convolutional encoder-decoder architecture for image segmentation. IEEE transactions on pattern analysis and machine intelligence. 2017. V. 39. № 12. P. 2481–2495. https://doi.org/10.1109/TPAMI.2016.2644615

  4. Benesty J., Chen J., Huang Y., Cohen I. Noise reduction in speech processing. Berlin, Heidelberg. Springer, 2009. V. 2. 230 p.

  5. Borovikov P.I., Antonov E.N., Dunaev A.G., Krotova L.I., Sviridov A.P., Fatkhudinov T.Kh., Popov V.K. Model of the destruction of bioresorbable polymers in aqueous media. Inorg. Mater. Appl. Res. 2018. V. 9. P. 649–654. https://doi.org/10.1134/S2075113318040056

  6. Brudfors M., Balbastre Y., Ashburner J., Rees G., Nachev P., Ourselin S., Cardoso M.J. An MRF-UNet Product of Experts for Image Segmentation. Proceedings of Machine Learning Research. 2021. arXiv preprint arXiv:2104.05495.

  7. Buzug T.M. Computed tomography: From photon statistics to modern cone-beam CT. Springer Science &Business Media. Berlin. 2008. 522 p.

  8. Chukalina M.V., Khafizov A.V., Kokhan V.V., Buzmakov A.V., Senin R.A., Uvarov V.I., Grigoriev M.V. Algorithm for post-processing of tomography images to calculate the dimension-geometric features of porous structures. Computer Optics. 2021. V. 45. № 1. P. 110–121. https://doi.org/10.18287/2412-6179-CO-781

  9. Gonzalez R., Woods R. Digital Image Processing. Up-per Saddle River, Prentice Hall. 2002. 1192 p.

  10. Gostick J.T., Khan Z.A., Tranter T.G., Kok M.D., Agnaou M., Sadeghi M., Jervis R. PoreSpy: A Python Toolkit for Quantitative Analysis of Porous Media Images. Journal of Open Source Software. 2019. V. 4. № 37. P. 1296–1299. https://doi.org/10.21105/joss.01296

  11. Jaccard N., Szita N., Lewis D. G. Segmentation of phase contrast microscopy images based on multi-scale local Basic Image Features histograms. Computer Methods in Biomechanics and Biomedical Engineering: Imaging & Visualization. 2017. V. 5. № 5. P. 359–367. https://doi.org/10.1080/21681163.2015.1016243

  12. Kurita T., Otsu N., Abdelmalek N. Maximum likelihood thresholding based onpopulation mixture models. Pattern Recognition. 1992. V. 25. № 10. P. 1231–1240. https://doi.org/10.1016/0031-3203(92)90024-D

  13. Myttenaere A., Golden B., Le Grand B., Rossi F. Mean Absolute Percentage Error for regression models. Neurocomputing. 2016. V. 192. P. 38–48. https://doi.org/10.1016/j.neucom.2015.12.114

  14. Ronneberger O., Fischer P., Brox T. U-net: Convolutional networks for biomedical image segmentation. Medical Image Computing and Computer-Assisted Intervention (MICCAI). Springer, LNCS. 2015. P. 234–241. https://doi.org/10.1007/978-3-319-24574-4_28

  15. Shah S.M., Gray F., Crawshaw J.P., Boek E.S. Micro-computed tomography pore-scale study of flow in porous media: effect of voxel resolution. Advances in Water Resources. 2016. V. 95. P. 276–287. https://doi.org/10.1016/j.advwatres.2015.07.012

  16. Shirani M.R., Safi-Esfahani F. Dynamic scheduling of tasks in cloud computing applying dragonfly algorithm, biogeography-based optimization algorithm and Mexican hat wavelet. The Journal of Supercomputing volume. 2021. V. 77. P. 1214–1272. https://doi.org/10.1007/s11227-020-03317-8

  17. Uvarov V.I., Borovinskaya I.P., Lukin E.S., Tsodikov M.V., Golubev K.B. Catalytically active cermet membrane for converting byproducts from the production of combustibles. Glass Ceram. 2014. V. 71. P. 270–274. https://doi.org/10.1007/s10717-014-9667-1

  18. Usanov M.S., Kulberg N.S., Morozov S.P. Experience of application of adaptive homophobic filters for computer tomograms processing. Journal of Information Technologies and Computing Systems. 2017. V. 2. P. 33–42.

Дополнительные материалы отсутствуют.