Изменения

Поиск повторов в ДНК на основе ОСАМ

4884 байта убрано, 23:24, 27 августа 2009
Нет описания правки
ОСАМ. Мысль простая: разложить сигнал по какому-нибудь классическому ортогональному базису, получить краткое описание, к тому же обладающее различными приятными свойствами. Обработать на основе описания сигнала. Применять можно в широком спектре задач распознавания. Свойства описания — «более важная» информация в первых коэффициентах; отсекая хвост, можно получать приближения сигнала; норма сохраняется; для неточных разложений есть мера точности разложения; и т. п. Т. е. есть хороший, проработанный, мат. аппарат.
Идея: применить ОСАМ к поиску повторов в ДНК, таким образом ускорив его. Как?! Во-первых, построить профиль последовательности, т.&nbsp;е. перевести её в длинный числовой вектор, выбрав w — окно профиля, и принимая за каждый элемент последовательности (количество пуринов в w-окрестности элемента) минус (количество пиримидинов в w-окрестности элемента). Далее, выбирая по N значений из полученной последовательности — 0..N-1, s..N+s-1, 2s..N+2s-1, … (s — шаг аппроксимации) и раскладывая получаемые вектора из N чисел по k коэффициентам некоторого базиса, получить «индекс» последовательности. k << N, потому и «индекс». Далее пробежаться по всем полученным описаниям (по индексу) обеих последовательностей (или одной и той же последовательности) и сравнить попарно все пары описаний (на похожесть). А что такое похожесть? Критериев похожести можно выработать массу, среди них можно найти устойчивые к масштабу и т.&nbsp;п., однако у нас всё довольно просто: <mathm>\frac{|a-b|}{|a|+|b|}</mathm>, где <mathm>|x|=\sqrt{\sum {x}_{i}^{2}}</mathm>. Такое вот «нормированное L<sub>2</sub>-расстояние». Здесь, кстати, можно выиграть от т.&nbsp;н. «принципа дискриминантности», который гласит очевидную вещь: что если <mathm>\frac{\sqrt{{\sum }_{i=0}^{k}{({a}_{i}-{b}_{i})}^{2}}}{|a|+|b|}> \mathrm{eps}</mathm> уже при k < n, то суммирование можно не продолжать, т.&nbsp;к. меньше сумма квадратов уже не станет. Итак, что мы получим от этого сравнения? Мы получим приближённые «близости» участков ДНК. Крупных или мелких, более или менее точное сравнение — это уже как захотим — для этого можно варьировать параметры. Задаём порог, можем пробежаться по результатам и сразу выявить «подозрительные на повторы» участки. Это есть важно, т.&nbsp;к. больше не нужно всё время искать повторы ВЕЗДЕ: сначала достаточно выявить крупные относительно похожие участки, а потом можно «увеличить масштаб» и выявить (или не выявить) точные координаты повторов. Кстати, единственное, для чего подход почти не подходит — для выявления «абсолютно точных» координат повторов. Это уже в «подозрительных» областях можно делать стандартными методами. Например, diffоподобным алгоритмом. :-)
Кстати, нужно использовать все современные возможности процессоров. Иначе будет обидно, если такую же программу написать на MATLAB’е и она — опа! — окажется быстрее в 5 раз. То есть нужно не забывать о многопоточности, не забывать об SIMD инструкциях, не забывать об аппаратном ускорении математических функций. Засчёт этого всего выигрываем в скорости ещё больше, реальная разница — в 10-20 раз (Core 2 Duo). Как?! Для многопоточности — голые нити (треды), никаких OpenMP! Так как это костылистая штуковина, приводит либо к сильному ухудшению структуры кода (причём фактическая логика получается аналогична голым тредам), либо к большим накладным расходам на распараллеливание — 5-15 %. Так что треды. Плюс библиотека Intel Integrated Performance Primitives для SIMD и аппаратного ускорения инструкций. А что это — IPP? А это такой векторный ассемблер, только на C. Библиотека, содержащая в себе оптимальные реализации большого спектра векторных операций (есть почти всё, что душе угодно — от сложений, умножений, корней и синусов, до узкоспециализированных функций ускорения декодирования аудио и видео, распознавания речи и т.д и т.п) для процессоров, имеющих различные расширения типа MMX / SSE1/2/3/4/5/+и т. д. Выражения над векторами там писать, к несчастью, нельзя, потому и получается код типа:
ippsCopy_64f(xn, wn, n);
Вот где-то примерно это всё и было реализовано. Есть относительно простая программа, есть относительно хорошая библиотека для абстрагирования от деталей реализации конкретных базисов, есть сами базисы — Чебышева 1 и 2 рода, Якоби, Лежандра, Лагерра, Эрмита, Фурье, ДКП, ДСП. Она работает и рисует красивые картинки. [показать пару картинок и закончить]. Кстати, по поводу того, а какой же базис лучше? Вообще они все дают очень похожие результаты… Пока что «лучше» всех Чебышев 1-го рода. А что вообще такое «лучше»? «Лучше» — чисто умозрительно это «больше соотношение сигнал/шум» (в результатах). Как измерить? Ну, например, при одинаковых параметрах окон и глубине разложения подобрать eps такое, чтобы общее количество «похожих» участков было примерно равно, и посчитать, например, среднюю длину повторов. Можно и медиану тоже. Чем больше, тем лучше — мы ведь хотим найти как можно более длинные повторы. Начинали реализовывать с Чебышева 1-го рода, потом пробовали Лежандра, потом думали, что Чебышев 2-го рода произведёт революцию и всё будет гораздо лучше, т.к весовая функция выпуклая, центр отрезка учитывается сильнее, края меньше. Революции не произошло, результаты сильно похожие на 1-го рода, местами получше, местами похуже. Формально — похуже. Дальше есть табличка с «попугаями» по разным базисам. Тестовые данные — часть генома мыши (не спрашивайте какая, я не знаю) длиной 1.5 млн нуклеотидов. Сравнение приводилось при примерно одинаковых количествах найденных участков, «подозрительных» на повтор — в районе 5000. При выбранных настройках минимальная длина участка, подозрительного на повтор — 3500 нуклеотидов. Какие выводы? Лидирует Чебышев 1 рода. Базисы ДКП, ДСП и Фурье дают до жути похожие на него, практически идентичные, результаты. С небольшим отставанием за ними следует Лежандр, за ним — Чебышев 2 рода, а базисы Эрмита и Лагерра не подходят для поиска повторов ''вообще — ''что есть логичный факт, т.&nbsp;к. они оба работают на бесконечном интервале либо (0, +бесконечность), либо от — до + бесконечности. Вариантов значения медианной длины было всего 2: 3500 (минимально возможная) или 10000, она отражает, фактически, чистое количество шума — мелких отрезков, и гласит, что приемлемый уровень шума дают… Ясно кто.
{| style<tab sep="tab"border-spacing:0;"| style="border-top:0.002cm solid #000000;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;1"|>| style="border-top:0.002cm solid #000000;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Eps| style="border-top:0.002cm solid #000000;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Среднее| style="border:0.002cm solid #000000;padding:0.097cm;"| Медиана |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Чебышева 1 рода| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .025| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| '''3978'''| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| '''10000''' |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Чебышева 2 рода| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .0285| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| 3882| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| 3500 |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| ДКП| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .025| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| '''3978'''| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| '''10000''' |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| ДСП| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .021| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| 3975| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| '''10000''' |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Фурье| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .025| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| '''3978'''| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| '''10000''' |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Эрмита| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .0015| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| 3502| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| 3500 |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Лагерра| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .0063| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| 3505| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| 3500 |-| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| Лежандра| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| .0225| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:none;padding:0.097cm;"| 3966| style="border-top:none;border-bottom:0.002cm solid #000000;border-left:0.002cm solid #000000;border-right:0.002cm solid #000000;padding:0.097cm;"| '''10000''' |}</tab>
[[Категория:Учёба]]