Определение регуляторных сегментов в геномах методами теоретического анализа последовательностей нуклеотидов ДНК, страница 5

Сегментацией конкретнойпоследовательности будем называть набор границ {b1 ,..., bK } разбивающийтекст T на сегменты. Каждому сегменту Ti сопоставляется вектор отсчетов –число встреченных букв каждого типа n = (n1 ,.., n L ) .При поиске оптимальной сегментации каждому сегменту соответствуетзависящая от отсчетов весовая функция – мера однородности сегмента.Полная весовая функция текста, состоящего из ряда сегментов равна суммевесовых функций каждого из сегментов. Для построения меры однородностисегмента используются методы байесовской статистики.

Для каждогосостава σ = {θ1 , θ2 ,..., θ L} задается функция плотности вероятности составов p(σ),Lкоторая определена на симплексе Ψ = {σ :θ k ≥ 0; ∑ θ k = 1} и удовлетворяетk =1условию нормировки ∫ dσp(σ ) = 1. Теорема Байеса позволяет вычислитьΨпостериорную плотность распределения вероятностей p( σ S ) на сегментеS: p( σ S ) =( )P S σ p(σ )P (S ), где P(S ) = ∫ dσ P(S σ ) p(σ ) это нормировка, называемаямаргинальным правдоподобием, зависящая только от вектора отчетов n и неменяющаяся при перестановках букв в последовательности S. Здесь( ) ∏i =1θ inP Sσ =L– вероятность реализации данной последовательности приизвестном составе, а ρ (σ ) - некоторое априорное распределение напространстве составов.

В задаче сегментации выбор априорногораспределения неявно означает определенное предположение о составеполимера, в наших расчетах используется однородное (неинформативное)распределение.В качестве меры однородости сегмента используется значения егомаргинального правдоподобия. Для однородного априорного распределения24маргинальное правдоподобие может быть получено аналитически (Liu,Lawrence, 1999):P(S ) = P(n ) =(L − 1)!( N + L − 1)!n1!...n L !.Вычисление оптимальной сегментации проводится согласно алгоритму,близкому к алгоритму Витерби (Forney, 1973). Пусть данапоследовательность S = s1s2s3… sN длины N, где si ∈ Σ . Пусть для каждогосегмента S(a,b) = sa… sb, a ≤ b может быть вычислен вес W(a,b).

В нашемслучае примем этот вес за равный ln ( P ( S ( a, b ) ) ) .Каждая конкретная сегментация R , имеющая m сегментов,определяется набором границ R = {k0 = 0, k1 ,..., km−1 , km = N } , где граница ki стоитмежду символами sk и sk . Определим вес всей сегментации Rim(i +1)F ( R ) = ∑ W k j −1 + 1, k j . Для вычисления оптимальной сегментации R , котораяj =1*максимизирует функционал F ( R ) , используется рекуррентный алгоритм,описанный в (Roytberg, Finkelstein, 1993).

Обозначим R* ( k ) оптимальнуюсегментацию фрагмента последовательности S (1, k ) , 1 ≤ k ≤ N . СегментацияR* (1) – тривиальна. В случае, если известны оптимальные сегментации всехучастков R* (1) ,..., R* ( k − 1) , можно найти оптимальную сегментацию участка[]R* ( k ) используя рекуррентное выражение: F ( R * ( k )) = max F ( R * (i )) + W (i + 1, k ) ,i = 0 ,..., k −1c начальным условием F(R*(0)) = 0.Более стабильной характеристикой является статсумма на множествесегментаций последовательности. Статсумма вычисляется как сумма весоввсевозможных сегментаций Z (N ) = ∑ ... ∑ Π(q1,... qN −1 ) , где индикатор qk равенq1qN−1единице, при границе после позиции k и нулю в других случаях, асегментация, заданная вектором границ q = ( q1 ,...qN ) , имеет вес (вероятность)Π(q).25Статистическая сумма может использоваться для оценки вкладаналичия границы на позиции k в разнообразные сегментации.

Для этогонеобходимо вычислить статистические суммы сегментаций фрагментовпоследовательности, ZL(k) и ZR(k) находящихся справа и слева от даннойграницы. Для вычисления этих значений могут использоваться следующиеk −1реккурентные формулы: ZL( k ) = ∑ eW ( j + 1, k − 1)j=0NZL( j ) , ZR( k ) = ∑ eW ( k , j ) ZR( j ) сj=kсоответствующими граничными условиям.Статистический вес границы, расположенной после позиции k, можетбыть вычислен как Π( qk = 1) =Z L ( k )Z R (N − k )Z (N ).

Границы с более высокимстатистическим весом разделяют сегменты, которые сильнее различаются посвоему составу.Глава 6. ВЫДЕЛЕНИЕ ХАРАКТЕРНЫХ МОТИВОВ ИЗВЫБОРОК НУКЛЕОТИДНЫХ ПОСЛЕДОВАТЕЛЬНОСТЕЙ СПОМОЩЬЮ АЛГОРИТМА ГИББСОВЫХ ВЫБОРОКВ главах 3 и 4 было рассмотрено вычисление статистическойзначимости встречаемости мотивов в ДНК. Экспериментальные данные,такие как SELEX и футпринтинг с ДНКазой I, позволяют получить мотивыдля каждого регуляторного белка; для этого используются методыопределения мотивов (motif discovery). В данной главе предлагается новыйметод определения мотивов, учитывающий симметрию участков ДНК,непосредственно взаимодействующих с регуляторным белком. Исходнойточкой послужил алгоритм Gibbs sampling (Lawrence, 1993).В алгоритме SeSiMCMC (сокращенно от «Sequence Similarity MarkovChain Monte Carlo) сохранено представление мотива в виде МПВ.

Каждаяисходная последовательность сегментируется на «мотив» и «фон»,порожденные разными вероятностными моделями. Фон моделируется какоднородная последовательность независимых случайных испытаний.26Критерием оптимальности сегментации считается ее максимальнаявероятность, вычисленная в рамках Байесовского подхода, при условииизвестной последовательности (ДНК данных) и отсутствия априорнойинформации о составе ДНК и предпочтения определенных нуклеотидов вкаждой позиции мотива. Выровненная часть набора последовательностейопределяет мотив. Оставшиеся невыровненными части всехпоследовательностей считаются фоном.Для выровненной части, для каждой позиции выравнивания можновычислить позиционные числа встречаемости нуклеотидов ni ,α , из которыхможно оценить позиционные вероятности qi ,α появления нуклеотида α впозиции i, i = 1..m , где m - это длина выравнивания, как q ( i, α ) =ni ,α + pαэтом для фонового распределения принимается оценка f (α ) =gα + bα.

Здесь MK+BM +B. При- это число выровненных последовательностей, K - это число всехневыровненных (фоновых) позиций во всех входных данных, а gα - полноечисло нуклеотидов типа α в невыровненных позициях. Псевдокаунты bαвыбираются пропорциональными частотам нуклеотидов во входных данных,в то время как их сумма B = ∑ bα ~ N , где N - это число входныхαпоследовательностей.Если принимается, что мотив имеет структуру прямого повтора, тоМПВ оценивается по формуле q ( i,α ) =ni ,α + n m +1 i +int ,α 2 + 2 ⋅ bα2 ⋅(M + B), в то время как дляпалиндромов (обратно-комплементарных повторов) формула принимаетвид:q ( i, α ) =ni ,r + nm+1−i ,α + bα + bα2 ⋅(M + B), где m - длина мотива и α — это нуклеотид,комплементарный α .Алгоритм действует следующим образом.

Вначале формируетсявыравнивание, состоящее из случайно выбранных участков некоторой длины,27по одному на последовательность. Исходная длина либо выбирается случайнов некоторых пределах, либо задается как параметр программы. Затемзаданные последовательности по очереди просматриваются (в цикле). Накаждом шаге выбирается «текущая» последовательность. Далее, повыравниванию, включающему в себя все последовательности, кроме«текущей», строится МПВ. Сегменты последовательностей, не вошедшие ввыравнивание, считаются порожденными из фонового распределения. Послетого, как оценены МПВ и фоновое распределение, вычисляется вероятностьполучить текущую последовательность, при условии, что мотив расположен впозиции k , и для каждой его позиции принята модель, взятая изсоответствующей колоноки МПВ:k −1k + m−1L − m+1i =1i =ki=k +mP (T | [ k ] , q, f ) = ∏ f ( ri ) ∏ q ( i − k + 1, ri ) ⋅P ( R | [ 0]) =∏ f ( r ), k ≠ 0;iL − m+1∏ f ( r ),i =1iздесь ri — это i -ый нуклеотид в последовательности T , а[ k ] , k = 1..( L − m + 1) обозначают событие: «сайт начинается с позиции k », [0]соответствует случаю отсутствия сайта (нулевая позиция).

Априорнаявероятность P ([0]) определяется пользователем и обозначает вероятностьтого, что последовательность из входных данных является шумом и не несетникакой биологической информации. Все ненулевые позиции имеют равныеаприорные вероятности.Апостериорная вероятность того, что сайт начинается в позиции kравна P ([ k ] | R, q, f ) =P (T | [ k ] , q, f ) P ([ k ])P ( T | q, f )=P ( T | [ k ] , q, f ) P ( [ k ] )P ( T | q, f ). Из этогоапостериорного распределения разыгрывается новая (возможно нулевая)позиция сайта. После этого полученная позиция сайта считается заданной длятекущей последовательности во всех последующих шагах цикла. Далеевыбирается следующая текущая последовательность, и цикл повторяется.Процесс последовательных итераций продолжается до тех пор, пока цепь28Анализаргументов иконфигурационного файлаОценкапредполагаемогоизменения качестванабора паттерновРешение: изменениепринимается или непринимается, вовтором случае новоесостояние повторяетстароеЧтение иверификациявходныхданныхСоздание иинициализациявнутреннихструктурВыбор типаследующего шагаизмененияпаттернов,розыгрыш новогопаттернаУвеличениесчётчика числашагов, егосравнение сзаданнымчислом шаговНовое состояниезаносится встатистикуВывод результатовАнализ собранной завремя работы статистикиРисунок 2.

Блок-схема работы алгоритма SeSiMCMC.полученных МЛАБП не сойдется (то есть, изменения от шага к шагу нестанут малыми). Мотив определяется из полного набора текущих положенийсайтов, который на рис.2 называется «паттерном».29Кроме местоположения мотива определяется оптимальная длина мотива идлины промежутка, если разрешены мотивы, состоящие из двух частей. Длякаждой длины мотива длина спейсера принимается как минимальноезначение, для которого достигается локальный максимум ИСП.