Как разобрать обезьяньи кишки на составные части. Изучаем цветовую деконволюцию

Как многие помнят, я работаю в лаборатории, где мы работаем с живыми и не очень организмами. Науку двигаем, короче. Обычно вперед. Иногда в качестве образцов нам достаются мертвые обезьяны, ткани которых потом идут на экспериментальные задачи. Выглядит обычно это крайне жизнерадостно. Раздается звонок в 11 часов вечера, и тебе сообщают, что в питомнике обезьянка убилась. Почти не поврежденная, соседи только сердце съели. Вздыхаем, лезем в расписание рейсов и едем в аэропорт. На месте тебе выдают нужные запчасти убиенной и складывают в прозрачный контейнер с консервационным раствором. В аэропорт с этим тащиться уже нельзя, так как ограничен провоз жидкостей. Идем на ж/д вокзал на экспресс до Краснодара. Милые девушки на контроле как правило приобретают восхитительный салатовый оттенок при виде медленно кружащихся органов в нежно-розовом растворе.
В-общем, привезли, нарезали все, что нужно ломтиками, покрасили… Но тут оказывается, что полученные исходники нужно обработать и посчитать в автоматическом режиме… Сразу хочу уточнить, что я врач-исследователь, а не профессиональный программист или математик. Поэтому, если что-то покажется ошибочным — буду рад правкам.

Формулируем задачу

Итак, в наличии есть серия срезов окрашенных двумя красителями. Голубоватый — гематоксилин Майера, коричневый — DAB хромаген. На этой иллюстрации пищевод обезьяны, а коричневым окрашены маркеры на панкератин. В принципе, это неважно. Главное, что у нас есть два красителя, которые могут одновременно окрашивать одни и те же структуры. Красить нужно одновременно обоими, иначе непонятно строение среза. Беда в том, что часто для подсчета объектов или определения площади структур с положительной реакцией на антитело нам нужно иметь изображение, где присутствует только один из пигментов. То есть задача может быть сформулирована примерно как "разделить обратно смесь кетчупа с майонезом". И здесь в игру вступает алгоритм цветовой деконволюции.

Цветовая деконволюция

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

I — значение интенсивности для данной длины волны. с — концентрация красителя, бета — коэффициент пропорциональности, который зависит от физических свойств данного вещества. То есть, по факту, нам нужно оценить тот вклад который вносит каждый из красителей в красный, зеленый и синий каналы пикселя. Для нашей задачи мы использовали алгоритм A.C. Ruifrok в модификации G. Landini, который написал плагин для открытого ПО ImageJ. ImageJ — софт для пакетной обработки изображения, написанный на Java. Крайне гибкий, позволяет обмазывать его самыми различными плагинами и писать под него скрипты. Сами мы ленивые, пусть за нас машина тонны изображений перемалывает. Но для начала нужно дать ей образцы, относительно которых нужно производить разделение.

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

if (myStain.equals("Custom DAB")){
		// This is the Custom DAB
			MODx[0]=0.66504073;
			MODy[0]=0.61772484;
			MODz[0]=0.41968665;

			MODx[1]=0.4100872;
			MODy[1]=0.5751321;
			MODz[1]=0.70785;

			MODx[2]=0.6241389;
			MODy[2]=0.53632;
			MODz[2]=0.56816506;

Теоретически можно взять готовые предустановленные лабораторные образцы. На практике получается, что в каждой лаборатории красители немного отличаются по оттенку, что делает крайне желательной калибровку на каждой партии. Для того, чтобы этот профиль появлялся в списке вариантов, придется засунуть это в исходный java-файл. Инструкция на странице разработчика. Чем лучше определены векторы красителей, тем ближе комплементарное изображение к белому — то, что остается после выковыривания наших пигментов. Итак, барабанная дробь и результат темной магии:

Пакетная обработка

Теперь нам нужно, чтобы все это безобразие работало в пакетном варианте. Мы же не хотим тыкать вручную в каждый файл? ImageJ позволяет писать пользовательские скрипты, что удобно для подобных задач. Собственно, создается обычный файл в формате txt, в котором описаны манипуляции с вашим исходником. У меня получилось примерно так:

dir1 = getDirectory("Choose Source Directory ");
dir2 = getDirectory("Choose Destination Directory ");
list = getFileList(dir1);
setBatchMode(true);
for (i=0; i<list.length; i++) {
showProgress(i+1, list.length);
open(dir1+list[i]);
imgName=getTitle(); 
run("Colour Deconvolution", "vectors=[Custom DAB]");
//close windows we don't need 
selectWindow(imgName + "-(Colour_3)"); 
close(); 
selectWindow(imgName +"-(Colour_1)"); 
title = getTitle(); 
saveAs("JPG", dir2 + title);
close(); 
selectWindow(imgName + "-(Colour_2)"); 
title = getTitle(); 
saveAs("JPG", dir2 + title);
close(); 
}

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

UPD. Автор алгоритма приводит достаточно интересные варианты применения — выявление надписей, залитых чернилами другого цвета (Beinecke Rare Book and Manuscript Library, Yale University)







UPD2. Вот эти милые бабуинчики очень любят убивать друг друга нахрен в иерархических стычках.