«Топологические особенности РНК-подобных молекул со случайной первичной структурой ...»

ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ БЮДЖЕТНОЕ УЧРЕЖДЕНИЕ НАУКИ

ИНСТИТУТ ХИМИЧЕСКОЙ ФИЗИКИ ИМ. Н.Н.СЕМЕНОВА

РОССИЙСКОЙ АКАДЕМИИ НАУК

На правах рукописи

ВАЛЬБА ОЛЬГА ВЛАДИМИРОВНА

Топологические особенности РНК-подобных молекул со случайной

первичной структурой Специальность 01.04.17 — Химическая физика, горение и взрыв, физика экстремальных состояний вещества Диссертация на соискание учёной степени кандидата физико-математических наук

Научный руководитель:

д.ф.-м.н., Аветисов В.А.

Москва – Оглавление Введение 1 Обзор литературы 1.1 Особенности пространственной структуры молекул РНК......... 1.2 Методы предсказания структуры РНК.................... 1.3 Случайная первичная структура РНК.................... 1.4 Термодинамические свойства........................ 1.5 Описание РНК структур случайными матрицами............. 2 Алгоритмы вычисления свободной энергии РНК-подобных структур 2.1 Выравнивание последовательностей..................... 2.2 Комплементарное связывание биополимеров................ 2.3 Связывание РНК с внутрипетлевым взаимодействием........... 2.4 Алгоритмы восстановления структуры................... 3 Свойства РНК структур со случайной последовательностью звеньев 3.1 Свободная энергия основного состояния.................. 3.2 Распределение длин петель в РНК-подобных структурах......... 4 Топология РНК-подобных молекул в зависимости от алфавита случайной первичной структуры 4.1 Зависимость свободной энергии РНК-подобных структур от алфавита. 4.2 Топологический переход в модели Бернулли................ 4.3 Аналитическая оценка критической точки топологического перехода в модели Бернулли................................ 4.3.1 Метод среднего поля......................... 4.3.2 Комбинаторная оценка........................ 4.3.3 Матричный подход.......................... 4.4 Переход случайной РНК в замороженное состояние, ограниченный топологическим переходом........................... 4.5 Другие модели нецелого алфавита...................... 4.5.1 Метод концентраций......................... 4.5.2 Коррелированная случайная последовательность......... 4.5.3 Рациональный алфавит........................ 5 Описание РНК-подобной структуры в терминах оптимизационной транспортной задачи 5.1 Оптимизационная транспортная задача................... 5.2 Модель случайных интервалов первичной структуры РНК-подобной молекулы..................................... 5.3 Топологические свойства РНК-подобных структур в модели случайных интервалов................................... 5.3.1 Численное моделирование...................... 5.3.2 Аналитическое описание....................... Заключение Список сокращений и условных обозначений Литература Введение Актуальность темы исследования. Структура важнейших биологических макромолекул, таких как дезоксирибонуклеиновые кислоты (ДНК), рибонуклеиновые кислоты (РНК) и белки, играет ключевую роль в их правильном функционировании в клетке.

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

Данная работа посвящена исследованию топологических свойств вторичной структуры молекул РНК-типа. Известно, что биомакромолекулы являются «слабо отредактированными случайными гетерополимерами» [1,2]. Более того, для ряда свойств распределение мономерных звеньев в первичной структуре, например, функциональных РНК можно считать случайным [3, 4]. В этом случае, модель случайной первичной структуры является базовой моделью, описывающей основной (нулевой) вклад в наблюдаемые физические явления. Основное внимание при этом сфокусировано на нетривиальной вторичной структуре РНК-подобных полимеров, для описания которой привлекаются разнообразные техники, в том числе, техники квантовой теории поля и моделей Изинга [5].

Цель работы заключается в описании топологических особенностей РНК-подобных последовательностей методами статистической физики и теории случайных процессов.

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

1. Разработать алгоритм вычисления свободной энергии РНК-подобной молекулы;

2. Исследовать статистические свойства распределения свободной энергии в ансамбле РНК-подобных структур со случайной последовательностью звеньев;

3. Изучить зависимость топологических свойств РНК-подобных структур от количества типов мономерных звеньев (алфавита), используемого в случайных первичных структурах.

Научная новизна работы заключается в следующем.

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

2. Впервые теоретически обнаружено критическое изменение топологии РНКподобных структур при переходе от двухбуквенного алфавита к трехбуквенному и проведена аналитическая оценка точки перехода в рамках комбинаторного и матричного описания.

3. Впервые установлена взаимосвязь между наблюдаемым критическим изменением топологии РНК-подобных структур и переходом в замороженное состояние, который обсуждался ранее в работах Т. Хва и Р. Бундшу.

4. Впервые показано, что описание топологии РНК-подобной структуры может быть сведено к оптимизационной транспортной задаче.

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

Методы исследования. В работе использовалось компьютерное моделирование, включающее вычисление свободной энергии основного состояния РНК-подобных молекул и предсказание соответствующих вторичных структур. В аналитическом рассмотрении широко использовалась теория случайных процессов, а также описание вторичной структуры РНК случайными матрицами.

Основные положения, выносимые на защиту:

1. алгоритмы описания вторичной структуры РНК-подобной молекулы и вычисления свободной энергии основного состояния, учитывающие внутрипетлевое взаимодействие;

2. свойства распределения свободной энергии в ансамбле РНК-подобных структур со случайной последовательностью мономерных звеньев;

3. зависимость топологических свойств РНК-подобных структур от используемого в первичной структуре числа различных мономерных звеньев (алфавита). Критическое изменение топологии РНК-подобных структур при переходе от двухбуквенного алфавита к трехбуквенному;

4. топологические свойства РНК-подобных структур с выбранным распределением расстояний между мономерными звеньями и потенциалом взаимодействия между мономерами, заданным выпуклой вниз функцией от расстояния.

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

Апробация работы. Основные результаты работы докладывались и обсуждались на 12 конференциях и 11 семинарах. По результатам диссертации опубликованы 6 статей в 5 ведущих российских и международных журналах.

Личное вклад автора заключается в развитии методов описания РНК-подобных молекул со случайной первичной структурой. Им были разработаны соответствующие алгоритмы вычисления свободной энергии РНК-подобных молекул. Все приведенные в работе расчеты и обобщение полученных результатов были выполнены автором лично.

Диссертация состоит из пяти глав и заключения. Первая глава содержит обзор литературных данных. Рассматриваются топологические особенности молекул РНК, приводятся существующие подходы к предсказанию вторичной структуры РНК. Отдельно обсуждаются термодинамические свойства РНК со случайной первичной структурой и матричный подход описания вторичной структуры РНК.

Во Второй главе приводятся алгоритмы описания РНК-подобной структуры и вычисления свободной энергии ее основного состояния. Формулируется вспомогательная статистическая модель, описывающая взаимодействия мономерных звеньев в РНКподобной структуре с петлевыми участками. Далее, в предложенной модели учитывается вклад внутрипетлевого взаимодействия мономеров и приводится соответствующий алгоритм динамического программирования для вычисления энергии такой иерархической структуры.

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

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

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

В Заключении представлены основные результаты диссертационной работы.

Глава Обзор литературы Глава посвящена обзору литературных данных. Обсуждаются особенности пространственной структуры молекул РНК, приводятся известные алгоритмы предсказания таких структур. Отдельно обсуждаются свойства РНК-подобных молекул со случайной последовательностью звеньев.

1.1 Особенности пространственной структуры молекул РНК — одна из трёх основных макромолекул, которые содержатся в клетках всех живых организмов. Так же, как ДНК, РНК состоит из длинной цепи, в которой каждое звено называется нуклеотидом. Последовательность нуклеотидов в цепи составляет первичную структуру РНК. Каждый нуклеотид состоит из азотистого основания, сахара (рибозы) и фосфатной группы. Последовательность нуклеотидов позволяет РНК кодировать генетическую информацию. Все клеточные организмы используют РНК для программирования синтеза белков, такая РНК называется матричной (мРНК).

Образование водородных связей между нуклеотидами обеспечивает вторичную структуру РНК. Азотистые основания в составе РНК могут образовывать водородные связи между цитозином C и гуанином G, аденином A и урацилом U. Такие пары называют комплементарными (Рис. 1.1) и впервые были обнаружены Дж. Уотсоном и Ф. Криком [6]. Помимо комплементарных пар, водородные связи могут образовываться между основаниями U–G (Рис. 1.1). Такие пары называются неканоническими (Wobble base pairs). Неканонические пары U–G влияют на вторичную и третичную структуры РНК и ее функции. В частности, было показано, что такие пары играют существенную роль в процессе кодон-антикодон связывания [7]. Теоретические расчеты показали, что энергия пары U–G сравнима с энергией основных пар C–G и A–U [8] (Табл. 1.1), однако геометрия пары отличается от канонических пар [9]. Гликозидный угол — угол, который образует связь между азотом N и сахарным остатком с плоскостью цикла сахара, одинаков для всех оснований в комплементарных парах и отличается от соответствующих углов в неканонической паре (Рис. 1.1).

Рис. 1.1 Отличие Уотсон–Криковских пар Таблица 1.1 Теоретическая оценка от неканонической пары G–U. свободной энергии и длины связи в Характерной особенностью структуры молекул РНК является то, что система комплементарных связей вторичной структуры представляет собой «клеверный лист»

(Рис. 1.2 (a)). Такая кактусообразная структура схематически может быть представлена набором вложенных дуг (Рис. 1.2 (в)), где дуга — связь между комплементарно связанными нуклеотидами в цепи. Псевдоузел — элемент вторичной структуры (Рис. 1.2 (б)) — образуется довольно редко и свойствен, в основном, длинным молекулам РНК. Псевдоузел соответствует пересечению дуг в арочном представлении вторичной структуры (Рис. 1.2 (г)). Псевдоузлы часто несут важную функциональную роль, например, было обнаружено что структура псевдоузла в теломерной РНК существенна для активности теломеразы [10].

Рис. 1.2 Клеверная структура РНК (a) и псевдоузел (б); (в) и (г) — арочное Задача предсказания вторичной структуры РНК достаточно сложна. Вторичная структура может содержать различные элементы, отличающиеся как по энергетическому вкладу в общую энергию биополимера, так и по энтропийному вкладу (Рис. 1.3) [11–14].

Рис. 1.3 (а) — 3D структура транспортной РНК дрожжей, полученная рентгеноструктурным анализом [15]; (б) — элементы клеверной структуры РНК [16].

1.2 Методы предсказания структуры РНК Наибольшую популярность приобрели методы предсказания вторичной структуры РНК, основанные на минимизации свободной энергии [11–13, 17, 18]. Основоположниками данного метода можно назвать М. Зукера и П. Стиглера [11]. В основе подхода лежит идея о том, что «правильная» вторичная структура РНК должна быть термодинамически наиболее стабильной и, следовательно, обладать наименьшей свободной энергией. При решении задачи минимизации энергии необходимы правила подсчета энергии для любой структуры и эффективный алгоритм минимизации энергии. На основе разнообразных экспериментальных данных [13] сделано много попыток построения правил подсчета свободной энергии и созданы достаточно эффективные алгоритмы, основанные на динамическом программировании [19]. Основное уравнение на статистическую сумму вторичной структуры РНК (Рис. 1.2 (а)) записывается как:

где, описывает статистический вес участка цепи с по мономер, а, определяется больцмановским весом контакта между и мономерами. Основное состояние определяется как:, = ln,, где и — константа Больцмана и абсолютная температура, соответственно. Так как энергия комплементарной связи превышает в десятки раз при комнатной температуре (Табл. 1.1), очень часто используют так называемое приближение нулевой температуры. В таком приближении, основное состояние определяется энергией взаимодействующих мономеров, тогда как, энтропией цепи можно пренебречь. Отметим, что выражение (1.1) может быть дополнено различными факторами, такими как минимальная длина петли, энергия стэкинга, различная энтропия структурных элементов РНК (Рис. 1.3 (б)). Особым случаем является предсказание пседвоузлов [20, 21], для которых разрабатываются отдельные алгоритмы с использованием динамического программирования. Методы, основанные на минимизации энергии, на сегодняшний день — наиболее часто используемые. Но, к сожалению, эти алгоритмы не являются надежными, и их точность сильно падает при увеличении длины последовательности. Также следует отметить, что в настоящее время еще не разработан подход, количественно оценивающий вероятность ошибочного предсказания РНК структуры.

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

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

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

Таким образом, предсказание вторичной структуры молекулы РНК по ее первичной — все еще открытый вопрос и исследования в этой области продолжаются [24–27].

Особое место среди таких задач занимают задачи о связывании РНК с биополимерами (белки, ДНК, РНК). Роль таких биополимеров как ДНК и РНК в механизмах клеточной регуляции общеизвестна. Их взаимодействие является одним из необходимых этапов клеточного цикла, связанного с хранением и передачей генетической информации. Помимо общеизвестных механизмов трансляции и транскрипции информации, основанных на ДНК–РНК связывании, исключительно важную роль играют РНК–РНК взаимодействия. Эти взаимодействия имеют ключевое значение для регуляции экспрессии генов [28, 29]. Молекулы РНК, посредством образования комплементарных пар, связываются с матричной РНК или ее участком и, тем самым, останавливают трансляцию генов с данной мРНК [28]. Молекулы РНК, участвующие в процессах данного типа, называются некодирующими РНК (нкРНК). Это название обусловлено тем, что они сами не транслируются в белки [29] и, следовательно, исключены из непосредственного процесса транскрипции.

Важная биологическая роль РНК–РНК взаимодействий обуславливает необходимость построения эффективного алгоритма, который бы позволил по первичным структурам молекул РНК, теоретически вычислять энергию связывания, а также предсказывать вторичную структуру такого комплекса. Эта задача тесно связана с проблемой выравнивания (alignment) двух произвольных линейных последовательностей типа ДНК.

Существенным отличием задачи выравнивания молекул РНК от аналогичной задачи для ДНК является наличие нетривиальной вторичной структуры у молекул РНК (Рис. 1.3).

Существует ряд подходов к определению энергии РНК-РНК взаимодействия [30–34].

Однако все они применимы в своем, достаточно узком семействе РНК последовательностей и «хорошо» работают только на конкретных примерах. Проблемы определения энергии РНК-РНК связывания аналогичны проблемам, возникающим в задачах предсказания вторичной структуры РНК, и эффективность того или иного алгоритма зависит от выбора факторов, которыми можно и нельзя пренебречь.

Конечно, ограничения того или иного метода могут оказаться существенными для предсказания структуры конкретной молекулы РНК, что, в свою очередь, может привести к неверным выводам о ее функции. Однако, для исследования статистических свойств случайных последовательностей РНК, т.е. цепочек со случайной первичной структурой, достаточно учесть основополагающие свойства полимера, — для РНК, это, в первую очередь, иерархическая вторичная структура типа клеверного листа, образующаяся согласно комплементарности азотистых оснований, и пренебречь теми, которые влияют, в большей степени, на структуру (и функцию) конкретной молекулы — псевдоузлами, минимальной длиной петли, стэкинг-взаимодействием.

1.3 Случайная первичная структура РНК Данная работа посвящена исследованию последовательностей со случайной первичной структурой. Известно, что биомакромолекулы являются «слабо отредактированными случайными гетерополимерами» [1, 2]. Более того, для ряда задач распределение мономерных звеньев в первичной структуре, например, функциональных РНК можно считать случайным [3,4]. Модель случайной первичной структуры является базовой моделью, описывающей основной (нулевой) вклад в наблюдаемые физические явления.

В работе рассматриваются статистические особенности вторичных структур длинных (от 1000 мономеров) случайных РНК последовательностей. Такие исследования играют важную роль, например, в понимании того, насколько «близки» или «далеки»

случайные РНК от реальных [35], какие свойства биополимера наиболее существенны для выполнения им определенной функции и, в конце концов, могли ли возникнуть функциональные РНК из случайных в ходе эволюции [36]. Случайные РНК представляют также довольно «удобную» систему для изучения термодинамических свойств реальных молекул. Исследование фазовых переходов [37–39], ответа цепочки на внешнюю силу [40, 41] основаны на модели случайной первичной структуры биополимера.

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

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

1.4 Термодинамические свойства В рамках модели случайной РНК последовательности были сделаны важные шаги в исследовании термодинамических свойств молекул РНК. Данные исследования важны не только для предсказания структуры и функции биополимеров, но также обширно используются для разработки методов скрининга экспериментальных данных для выявления генетических маркеров заболевания [42], секвенирования одиночных нуклеотидных полиморфизмов, выбора оптимальных условий для экспериментов по гибридизации и клонированию [43, 44]. Кроме того, разработка ДНК-чипов для быстрого скрининга и секвенирования основана на способности предсказывать термодинамическую устойчивость комплексов, образованных олигонуклеотидными зондами [45, 46].

С пионерских работ Бундшу и Хва [37,38], несколько авторов занимались исследованиями термодинамических свойств случайных РНК [39,47–49]. К настоящему времени, принято считать, что в этой системе имеет место фазовый переход в «замороженное»

состояние при низких температурах. Основываясь на репличном анализе, Лассиг и Визе, [50] и Давид и Визе [51] сформулировали задачу о переходе в терминах теории поля. Ниже приводятся доводы Бундшу и Хва, доказывающие существование фазового перехода и обсуждаются характерные свойства разных фаз.

В зависимости от температуры, случайная РНК находится в одной из фаз: i) «расплавленная» высокотемпературная фаза (molten phase) или ii) «замороженная» низкотемпературная фаза (glass phase). В высокотемпературной фазе большую роль играет энтропия цепочки, нежели порядок мономеров в первичной структуре. Данная фаза хорошо описывается в модели гомополимера, комплементарное связывание не играет роли, и эффективно можно заменить все мономеры мономерами одного типа. Низкотемпературная фаза, наоборот, определяется, в первую очередь, первичной структурой цепочки, то есть, основной вклад в свободную энергию обусловлен комплементарным связыванием мономеров. Такую фазу принято характеризовать замороженным беспорядком [37, 38]. Температура, при которой РНК переходит из одной фазы в другую, называется температурой фазового перехода и в литературе обозначается.

Был предложен следующий подход к определению температуры фазового перехода. Рассмотрим пару мономеров, чье взаимодействие приводит к образованию петли наибольшего размера, т.е. нуклеотидов с номером 1 и /2 по цепи для последовательности длиной (Рис. 1.4(а)). Определим энергию выигрыша данного контакта, которая определяется как () = ln 1,/2, где 1,/2 — вероятность связывания 1 и /2 мономера цепи. Данную энергию называют энергией пинча, и из выражения для статистической суммы цепочки (1.1), легко видеть, что:

Вероятность образования контакта между мономерами в высокотемпературной фазе пинча, таким образом, линейно зависит от температуры. Температура, при которой нарушается линейная зависимость (), и есть температура фазового перехода. В численном моделировании температуру перехода обычно определяют следующим образом. Зависимость () от длины случайной последовательности РНК аппроксимируют прямой и строят зависимость угла наклона ( ). В высокотемпературной фазе с хорошой точностью ( ) = 3 (Рис. 1.5). В низкотемпературной фазе, в отличие от высокотемпеa) Рис. 1.4 Вычисление энергии пинча: разделение цепочки на две половинки ограничивает число возможных конфигураций (a). Полностью комплементарные участки, один из которых лежит между 1 и /2, а второй — между (/2 + 1) и (б) могут быть найдены почти для любой случайной РНК. Комплементарное связывание ограничивает возможные конфигурации (в). Молекула разделяется на две петли, в каждой из которых образование связей происходит независимо [38].

ратурной, выигрыш зависит, в первую очередь, от первичного беспорядка (структуры) последовательности (Рис. 1.4(б,в)). Разрыв контакта между 1 и /2 мономером определяется не столько энтропийным фактором — насколько близки или далеки данные мономеры, а скорее энергетическими, т.е. средней энергией на мономер, числом несвязанных мономеров в цепочке. Зависимость наклона ( ) = 2 нарушается (Рис. 1.5).

Низкотемпературная фаза характеризуется линейным ростом энергии пинча с уменьшением температуры. Точка пересечения двух прямых определяет температуру фазового перехода. В работе [48] было высказано предположение о том, что в низкотемпературной фазе, энергия пинча зависит от логарифма длины не линейным образом (Рис. 1.5), а квадратично.

Переход между расплавленной и замороженной фазой относится к непрерывным переходам второго рода [50]. Было показано, что температура перехода непосредРис. 1.5 Зависимость наклона ( ) уравнения (1.4) от температуры для случайной последовательности РНК [38]. Температура представлена в единицах энергии ственно связана со средним количеством несвязанных мономеров в основном состоянии [38]. Аналитическая оценка температуры перехода на порядок отличается от экспериментально полученной [38].

Высокотемпературная и низкотемпературная фазы обладают разными скейлинговыми свойствами. Одна из величин, которая представляет интерес, — характерный размер структуры РНК. Под характерным размером понимают высоту соответствующей диаграммы в арочном представлении (Рис. 1.6) Было показано численно [38] и затем подтверждено аналитически [51], что низкотемпературная фаза характеризуется степенной зависимостью от длины последовательности со степенью 0.64, что близко к 0 = 2/3, и указывает на класс универсальности Кардара-Паризи-Жанга [52], характерного для таких процессов как, например, рост поверхности и баллистическая депозиция [53]. В высокотемпературной фазе численный эксперимент дает степень 0.54 [38], что находится в согласии с ожидаемой Рис. 1.6 Характерный размер структуры РНК в арочном представлении (а). Размер структуры определяется количеством пар, которые нужно разбить в максимально ит отметить, что значительную роль в исследовании скейлинговых свойст случайной РНК сыграл матричный подход к описанию структуры. Построенная полевая теория перехода позволила также говорить о том, что переход происходит через образования зародышей в расплавленной фазе [50]. В следующем разделе подробно представлены основные положения описания структуры РНК случайными матрицами.

1.5 Описание РНК структур случайными матрицами Для простоты предположим, что случайный полимер обладает бесконечной гибкостью, таким образом, можно пренебречь стерическими ограничениями и говорить, что любые парные взаимодействия в цепочке возможны [54]. Статистическая сумма такой последовательности длиной в этом случае может быть представлена в виде:

где, = (, / ) обозначает статистический вес контакта (, ) с соответствующей энергией, ; < > обозначает все пары <, < > — четверки < < < и т.д. Суммирование ведется по все возможным контактам в цепи. Как было показано в [55], каждый член соответствует своей арочной диаграмме (Рис. 1.2). В этом представлении нуклеотиды — это точки, ориентированные на горизонтальной оси в направлении от 5‘ к 3‘ концу молекулы и каждая комплементарная пара — арка между взаимодействующими основаниями. Диаграммы, состоящие из непересекающихся арок называют планарными. Такие диаграммы соответствуют кактусообразным структурам РНК. Структуры типа псевдоузлов в диаграммном представлении, как уже упоминалось, соответствуют пересечению арок. Основная идея матричного анализа структур РНК заключается в следующем [55]. Рассмотрим интеграл по случайным матрицам Здесь, где пробегает значения от = 1 до =, обозначает -ую случайную эрмитову матрицу, размера и (1 + ) — упорядоченное произведение таких матриц: (1 + 1 )(1 + 2 )...(1 + ). Нормировочный множитель и — симметричная матрица с элементами,. Интеграл 1.7 можно оценить используя теорему Вика. В результате, для больших интеграл (1.7) можно представить, как Взаимосвязь полученной функции со статистической суммой (1.6) очевидна. При = функции совпадают, для > 1 слагаемые (1.9) содержат информацию о топологии конфигураций. Все планарные диаграммы описываются членом (1) (1.9), более высокие порядки разложения 1/ 2 соответствуют структурам РНК с псевдоузлами. Классификация псевдоузлов, возникающая из выражения (1.9) описана в [56]. В общем случае для конкретной последовательности РНК, описываемой матрицей с элементами,, задача определения всех возможных конфигураций оказывается очень сложной. Для точного описания топологии вторичных структур РНК пользуются рядом упрощений. Так, предположение, что все элементы матрицы, равны между собой, =, позволяет вывести общую формулу для статистической суммы, причем каждое слагаемое несет информацию и о топологии структуры и о количестве контактов в ней (Табл. 1.2). В этом случае, многомерный интеграл (1.7) может быть сведен преобразованиями ХаббардаСтратоновича к одномерному, включающему спектральную плотность гауссовой матрицы [57]. Так как выражение для спектральной плотности хорошо известно из теории случайных матриц [58], интеграл (1.7) можно вычислить точно. Так, для четырехбуквенной последовательности статистическая сумма 4 ( ) = 1 + 6 + 2 2 + 2 / 2, а соответствующие разрешенные конфигурации представлены на Рис. 1.7.

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

1/ 2 –разложение статистической суммы ( ) можно представить где, описывает количество конфигураций определенного типа ( = 0 — планарные диаграммы, = 1 — конфигурации с одним псевдоузлом) c арками (Рис. 1.7).

Разложение ( ) по степеням большим чем 1/ 2 описывает конфигурации со сложными псевдоузлами. Теория таких структур описана в [60].

Таблица 1.2 Разложение статистической суммы (1.7) по степеням 1/ 2 для разных Глава Алгоритмы вычисления свободной энергии РНК-подобных структур Данная глава посвящена разработанным алгоритмах описания вторичной структуры биополимеров и вычисления их свободной энергии.

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

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

2.1 Выравнивание последовательностей Задача о выравнивании двух последовательностей – это задача нахождения эффективного алгоритма поиска наибольшей общей подпоследовательности (НОП) двух произвольных линейных последовательностей. Данная проблема является одной из ключевых задач вычислительной эволюционной биологии. В частности, она позволяет судить о том насколько далеко (в эволюционном смысле) разошлись друг от друга два рассматриваемых гена и какие гены могут являться их общими предками [61, 62]. Задача об НОП широко исследовалась в биологии [63–65], компьютерных науках [66–69], теории вероятности [70–75] и позже в статистической физике [37, 76–78].

Задача о поиске НОП двух последовательностей формулируется следующим образом. Рассмотрим две произвольные последовательности (в качестве примера рассматриваются последовательности РНК, составленные из 4-х буквенного алфавита A, C, G, U):

1 = {A, C, G, C, U, A, C} длины = 7 и 2 = {C, U, G, A, C} длины = 5. Далее, везде под алфавитом подрузамевается количество различных мономерных хвеньев в первичной структуре. Общая подпоследовательность – это подпоследовательность, содержащая буквы (нуклеотиды) как первой, так и второй последовательности, причем подпоследовательность необязательно содержит буквы, идущие непосредственно друг за другом. Так, например, для двух последовательностей 1 и 2 можно выделить несколько различных общих подпоследовательностей, например, {C, U, A, C} или {G, A, C} – обе эти подпоследоавтельности содержатся в 1 и 2, и являются для них общими.

Число возможных общих подпоследовательностей с ростом длин и полимеров растет экспоненциально. Алгоритм для определения оптимального выравнивания двух последовательностей впервые был сформулирован в [79]. В наиболее общем смысле каждое выравнивание двух последовательностей характеризуется числом совпадающих и несовпадающих букв и числом пропусков (делеций) в выравненных последовательностях. Для каждого выравнивания можно ввести весовую функцию (cost function), имеющую значение энергии [79, 80]:

В формуле (2.1) match, mis и gap – число пар совпадающих букв, число пар несовпадающих букв и число делеций в рассматриваемом выравнивании, соответственно.

Величины и – это вклады в весовую функцию от пары несовпадающих букв и делеции; вклад от пары совпадающих нуклеотидов, без потери общности, можно считать равным 1. В таком представлении функция удовлетворяет очевидному закону сохранения:

Используя (2.2), формулу (2.1) можно переписать в виде:

где Здесь интерес представляет область 0 1, так как, случай < 0 неотличим от = 0, а случай > 1 соответствует тому, что «несовпадения» более выгодны, чем «совпадения» и может быть учтен простым переопределением этих понятий. Заметим, что, хотя предлагаемая теория применима ко всему доступному интервалу значений, все численные результаты настоящей работы получены для случая = 0, который представляется наиболее физически осмысленным. Задача поиска НОП заключается в определении выравнивания с максимальным значением весовой функции.

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

где Выражения (2.5)—(2.6) имеют следующий смысл. Начиная с левых концов последовательностей, на каждом шаге выбирается такое положение букв в выравнивании, которое вносит наибольший вклад в функцию. Члены в (2.5) соответствуют трем возможным ситуациям: пропуску буквы в первой последовательности, пропуску во второй последовательности и случаю, когда -ая буква первой последовательности выравнена с -ой буквой второй последовательности.

2.2 Комплементарное связывание биополимеров Цель работы заключается в разработке статистического алгоритма вычисления весовой функции, которая бы характеризовала «похожесть» двух заданных последовательностей со сложной вторичной структурой типа РНК. Эта функция должна включать как энергетический вклад от непосредственного взаимодействия мономеров друг с другом, так и энтропийный вклад, обусловленный наличием ансамбля пространственных конформаций макромолекул. При этом постараемся, по возможности, остаться в рамках статистической физики и избежать неконтролируемых эвристических соображений, апеллирующих к опыту, полученному в результате анализа экспериментальных данных.

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

Рассмотрим вспомогательную статистическую модель, описывающую взаимодействие двух линейных полимеров с произвольными первичными последовательностями.

Пусть длины этих последовательностей, измеренные в единицах мономерных звеньев, равны и, соответственно. Каждый мономер может быть выбран из различных мономеров,,,,... (Для последовательностей РНК = 4). Мономеры первой последовательности могут образовывать связи с мономерами второй последовательности.

В молекулах РНК такие связи образуются согласно комплементарности азотистых оснований (1.1). Будем считать энергию связи между комплементарными нуклеотидами равной, а энергию между некомплементарными равной, где и — некоторые положительные величины (|| > ||). Предположим также, что некоторые части полимеров могут образовывать петли. На Рис. 2.1 схематически представлено взаимодействие двухбуквенных полимеров. Очевидно, что петли соответствуют делециям в задаче о выравнивании двух последовательностей.

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

Рис. 2.1 Взаимодействие мономерных звеньев в РНК-подобной структуре с петлевыми участками как выравнивание соответствующих последовательностей При низких температурах, можно представить как:

Смысл данной формулы очевиден: начиная с левого конца последовательностей (Рис. 2.1), находим первый существующий контакт между -м мономером первой цепи и -м мономером второй, а далее суммируем по всем возможным расположениям этого контакта. Статистические веса связей, определяются энергией контакта между -ым и -ым мономерами:

Здесь и далее,. Легко проверить, что статистическая сумма вида (2.7) удовлетворяет рекуррентному соотношению:

В свою очередь, статистическая сумма связана со свободной энергией комплекса, и температурой известным соотношением, = exp{, / }. Будем интересоваться значением свободной энергии с точностью до знака, тогда для величины, =,, переходя в уравнении (2.9) к пределу 0, получим:

Формулу (2.10) можно переписать в виде:

где введено обозначение:

Принимая + за единицу энергии, перепишем формулу (2.11) в виде:

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

Отметим, что предложенная выше модель является лишь первым приближением к описанию комплексообразования биополимеров. Известно (см., например, [2]), что для точного количественного описания такого связывания, например, двойной спирали ДНК необходимо учесть еще ряд факторов. Во-первых, не учтены так называемые «петлевые факторы»: при образовании петли возможные конформации полимера ограничены условием, что ее концы обязаны сойтись в одной точке пространства, поэтому образование каждой петли приводит к снижению энтропии комплекса. Во-вторых, в реальной ДНК имеется выраженная кооперативность образования связей: вероятность образования связи выше, если соседние мономеры также образуют связь. В-третьих, не учтено, что гибкость полимера конечна и, тем самым, существует ограничение на минимальную длину петли. И наконец, не было принято во внимание то обстоятельство, что комплементарные пары и имеют различную энергию связи и, что помимо комплементарных пар, возможно образование неканонических пар (см. 1.1).

Обобщение выражений (2.7)–(2.14) с учетом кооперативности образования связи, минимальной длины петли и различной энергией комплементарных связей — задача вычислительно сложная, но не требующая качественного изменения предложенного формализма, т.к. эти факторы влияют только на локальные свойства полимерных цепей. С другой стороны, петлевой фактор — характеристика нелокальная, зависящая от расстояния между мономерами, образующими связи и в этом случае нельзя описать состояние комплекса уравнениями динамического программирования, вида (2.7).

Однако, поскольку петлевые факторы имеют энтропийную природу, в пределе низких температур ( 0) их вклад в свободную энергию гетерополимерного комплекса становится пренебрежимо мал. Ситуация усложняется, если сами петли могут образовывать вторичную структуру (т.е. если внутри петли имеет место взаимодействие между мономерами), а именно такая ситуация типична для последовательностей РНК.

В этом случае энергетический вклад от вторичной структуры петли сохраняется и в пределе нулевой температуры, и его учет становится необходим.

2.3 Связывание РНК с внутрипетлевым взаимодействием В этом разделе обобщается модель взаимодействия двух сополимеров на случай, когда возможно комплементарное связывание внутри петель комплекса. Будем рассматривать иерархические структуры петель типа клеверного листа (Рис. 1.2(a)), структуры типа псевдоузлов (Рис. 1.2(б)) в данной работе не рассматриваются. Как и в предыдуРис. 2.2 Диаграмма для вычисления статистического веса последовательности.

щем параграфе для простоты не будем учитывать кооперативность образования связей и различие в энергиях комплементарных пар. Однако, как уже указывалось, модель может быть обобщена с учетом этих факторов. Согласно [41] можно переписать выражение (2.7) для статистической суммы, двух взаимодействующих сополимеров в виде:

где, и, обозначены статистические веса участков (с -го нуклеотида до -го) первой и второй последовательности, соответственно, удовлетворяющие уравнениям [81]:

Эти уравнения отвечают за топологию кактусообразной структуры, свойственной молекулам РНК, диаграмма, описывающая такие структуры представлена на Рис.2.2. Коэффициенты, – это константы, описывающие взаимодействие внутри последовательности, аналогичные,. Суммирование по ведется от + 1 + до для того чтобы исключить петли длиной меньше мономеров. В последующих вычислениях, как правило, предполагается, что =0, также обсуждается случай = 3. Напомним еще раз, что так как интерес представляет низкие температуры, можно пренебречь вкладом, связанным с потерей энтропии при образовании петель.

Сложную систему уравнений на статистические веса петлевых участках,, = 1, 2 (2.16) можно решить следующим образом. Для каждой из последовательности РНК можно построить матрицу, (, )-й элемент которой определяет статистический вес участка, начинающейся с -го нуклеотида и заканчивающейся -м. Таким образом, статистические веса всех возможных петель описываются матрицами размера для первой последовательности и для второй. Из граничных условий (2.16) можно однозначно определить элементы,+1. Из (2.16) следует, что элементы последующих субдиагоналей,+ зависят только от элементов предыдущих субдиагоналей,+ матрицы:

Определенные таким образом матрицы статистических весов () всех возможных петель позволяют вычислить статистическую сумму взаимодействия двух РНК с внутрипетлевым взаимодействием (2.15).

Как и в случае связывания последовательностей без петлевых участков, можно выполнить переход к пределу нулевой температуры – см. выражения (2.7)-(2.14). Элементы матрицы свободной энергии при этом можно представить в виде:

где, = lim ln, ( = 1, 2) имеют смысл с точностью до знака свободных энергий петлевых участков последовательностей с -го нуклеотида по -й,, – (, )ый элемент суммы (2.15), который в пределе нулевой температуры есть:

Элемент, описывает энергию комплекса взаимодействующих РНК, не имеющих контакта правее пары (, ). Из (2.17) следует, что функции, удовлетворяют:

здесь величина, – величина, как в (2.14),, – аналогичная величина, описывающая взаимодействие внутри петель. На свободную энергию накладываются граничные условия, как это следует из (2.15):

Таким образом, для того, чтобы вычислить энергию основного состояния комплекса двух взаимодействующих РНК, необходимо построить матрицы (1) и (2) и, далее, применяя (2.18)-(2.19), определить элементы матрицы.

Отметим, что выражения (2.17), (2.20) можно использовать для непосредственного вычисления свободной энергии основного состояния одноцепочечной РНК.

2.4 Алгоритмы восстановления структуры В данном разделе показывается, как алгоритм для вычисления энергии основного состояния может быть применен для восстановления структуры. Рассмотрим сначала комплексообразование двух линейных сополимеров.

Нахождение НОП двух линейных последовательностей В отличие от (2.13) будем теперь интересоваться, не количеством мономеров в наибольшей общей подпоследовательности, а ее составом, т.е. определением, из каких мономеров состоит общая НОП. Отметим, что, вообще говоря, задача может иметь множество решений, т.е., основное состояние может быть вырождено. Алгоритм, который будет описан ниже, позволяет определить все возможные НОП двух цепочек. Здесь и далее, в численном моделировании использовались параметры || = 1 и || = 0. В таком рассмотрении функция (??),(??) совпадает с количеством комплементарных связей в структуре основного состояния.

Рассмотрим действие алгоритма на конкретном примере. Возьмем две последовательности:

(очевидно, что = = 6), и построим матрицу с, = 1, если –й мономер первой последовательности комплементарен –му мономеру второй последовательности и, = 0, в противном случае (см. Рис. 2.3(а)). Далее, строим матрицу, используя рекурсивный алгоритм (2.13)-(2.14) (см. Рис. 2.3(б)). Нижний правый элемент этой матрицы 6,6 = 4 соответствует количеству комплементарных связей в основном состоянии комплекса. Теперь, чтобы установить, какие конкретно мономеры образуют связь, нам нужно, восстановить каждый шаг алгоритма (2.13). Вся информация содержится в матрице. Действительно, сравним (, )–элемент матрицы с соседними ему элементами 1,1, 1,,,1 и:

1. если 1,1 = max [1,1, 1,,,1 ], тогда в оптимальном выравнивании -й нуклеотид первой последовательности связан с -м мономером второй;

2. если 1, = max [1,1, 1,,,1 ], то в оптимальном выравнивании -й нуклеотид первой цепочки не участвует в комплексообразовании и таким образом, является пропуском (петлевым мономером);

3. и, наконец, если,1 = max [1,1, 1,,,1 ], то пропуском является -й мономер второй последовательности.

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

Рис. 2.3 Алгоритм восстановления связей при взаимодействии двух РНК с петлевыми участками: матрица возможных контактов (а) и матрица (б), построенная согласно (2.13)-(2.14); оптимальные пути на матрице и соответствующие им Структура комплекса с внутрипетлевым взаимодействием Более сложной является процедура восстановления структуры комплекса с внутрипетлевым взаимодействием. Схематически алгоритм определения контактов в оптимальной конфигурации показан на Рис. 2.4.

Рис. 2.4 Алгоритм определения оптимальной конфигурации комплементарного связывания РНК-подобных молекул с внутрипетлевым взаимодействием.

Начальным элементом, как и в задаче без петлевых участков, выбирается нижний правый элемент матрицы,. Если, > 1, + 1, (см. (2.18)), то информация о контакте в оптимальной конфигурации содержится в матрице (2.19). Отметим, что каждая пара (, ) характеризуется своей матрицей. Максимальный элемент матрицы, говорит о контакте между -м мономером первой последовательности и -м нуклеотидом второй (Рис. 2.4). Далее описанная процедура повторяется. Как и в предыдущем параграфе, рассмотрим действие алгоритма на конкретном примере:

Для последовательностей 1 и 2, на Рис. 2.5 приведены соответствующие матрицы (2.20)-(2.20). Элемент 7,7 = 6 показывает, что в оптимальной конфигурации содержатся 6 комплементарных пар. Чтобы, установить какие именно мономеры образуют пары, согласно разработанному алгоритму, рассматриваем матрицы для последовательно устанавливаемых контактов. Так, максимальный элемент 7,7 для 7, показывает, что 7-ой мономер 1 и 7-ой мономер 2 образуют пару. На следующем шаге рассматриваем матрицу для элемента 6,6 (2.19). Если матрица содержит несколько одинаковых максимальных элементов, это свидетельствует о вырожденности основного состояния. В данном случае процедура повторяется для каждого из элементов. Отдельно восстанавливается структура комплементарных связей внутри петель комплекса. Отметим, что эта задача идентична восстановлению структуры отдельной цепочки РНК. Алгоритм восстановления структуры петли основан на (2.20) и заключается в установлении, какая пара мономеров (, ) обеспечивает наибольший вклад в выражение для. На Рис. 2.5 приведены две возможные конфигурации рассматриваемых последовательностей 1 и 2.

Рис. 2.5 Алгоритм определения оптимальной конфигурации связывания РНК с внутрипетлевым взаимодействием: матрицы контактов внутри последовательностей (а, б) и между ними (в); матрицы 1 (г) и 2 (д), и (е), вычисленные по (2.18)-(2.20);

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

Рис. 2.6 Комплементарное связывание двух РНК: с петлевыми участками (a), с внутрипетлевым взаимодействием и минимальной длиной петли = 0 (б), и = 3 (в).

Глава Свойства РНК структур со случайной последовательностью звеньев В данной главе обсуждаются свойства распределения свободной энергии основного состояния в ансамбле РНК-подобных молекул со случайной первичной структурой.

Также, приводятся результаты для распределения длин петель в РНК-подобных структурах и обсуждаются аналитические модели их описания.

3.1 Свободная энергия основного состояния Связывание двух РНК с петлевыми участками Задача поиска оптимальной конфигурации линейного выравнивания случайных последовательностей неоднократно рассматривалась в литературе (см., например, [82,83]) в рамках так называемой модели «бернуллиевского сравнения», т.е. в предположении о том, что матричные элементы, (2.14) являются независимыми случайными величинами, принимающими значения 1 с вероятностью = 1 и 0 с вероятностью = 1, где — алфавит, используемый в случайной первичной структуре полимера. В работе [83] было показано что для длин последовательностей, 1 распределение энергии основного состояния имеет вид:

где – случайная величина с распределением Трейси–Видома ( = 1.7711... и 2 2 = 0.8132...)(более подробное описание этого распределения можно найти, например, в обзоре [84]). При =, оптимальная конфигурация характеризуется:

Флуктуации свободной энергии подчиняются: [83]:

Показатель 1/3 является типичным для стохастической динамики сильно коррелированных систем и относится к классу универсальности Кардара-Паризи-Занга (Kardar– Parisi–Zhang (KPZ)) [52].

Результаты численного моделирования распределения свободной энергии основного состояния для ансамбля случайных первичных структур РНК представлены на Рис. 3.1.

Угловой коэффициент прямой 0.65 (Рис. 3.1(a)), что хорошо согласуется с веливычисленной по формуле (3.2). Для флуктуации энергии полученный наклон 0.34 (Рис.3.1(б)) также близок к значению 1. Таким образом, уравнение (3.2), полученное в приближении бернуллиевского сравнения, удовлетворительно описывает численно наблюдаемую зависимость энергии основного состояния при связывании сополимеров с петлевыми участками от длины случайных цепей.

Связывание двух РНК с внутрипетлевым взаимодействием Аналогичный анализ был проведен и для двух последовательностей, образующих структуру с внутрипетлевым взаимодействием и минимальной длиной петли = 0. Соответствующие графики зависимости свободной энергии и флуктуации энергии представлены на Рис. 3.2. Как и для взаимодействия с петлевыми участками,, () = при 1 (Рис. 3.2), но угловой коэффициент прямой 0.92 гораздо выше, что обусловлено взаимодействием нуклеотидов внутри петель. Зависимость флуктуации энергии основного состояния остается такой же (см. Рис. 3.2(б)).

Рис. 3.1 Взаимодействие РНК с петлевыми участками: зависимость среднего значения свободной энергии основного состояния, (а) и флуктуации энергии (б) от длины случайной последовательности. Усреднение проводилось по ансамблю из случайных пар последовательностей для каждого значения длины.

Рис. 3.2 Связывание РНК с внутрипетлевым взаимодействием: зависимость энергии основного состояния, (а) и флуктуации свободной энергии (б) от длины случайной последовательности. Усреднение проводилось по ансамблю из случайных пар последовательностей для каждого значения длины.

Рис. 3.3 Иерархическая модель связывания двух полимеров с внутрипетлевым взаимодействием. Петли первого ( = 1), второго ( = 2) и третьего ( = 3) Оценим аналитически величину коэффициента в зависимости свободной энергии от длины цепи для внутрипетлевого взаимодействия (Рис. 3.2). Будем рассматривать комплекс, который образуют две случайные последовательности РНК, как структуру, состоящую из петель различных иерархических уровней, занумерованных индексом (см. Рис. 3.3).

Каждую петлю -ого иерархического уровня можно рассматривать как комплекс двух взаимодействующих подпоследовательностей из которых она состоит. Из выражения (3.1) следует, что наибольший вклад в свободную энергию наблюдается для комплекса, состоящего из двух последовательностей равной длины, =. Это позволяет оценить сверху свободную энергию петли как свободную энергию двух взаимодействующих половинок этой петли. Представление комплекса двух молекул РНК в виде иерархической структуры позволяет использовать идеи ренормализационной группы [85]. А именно, комплексы -ого иерархического уровня содержат петли, которые будем считать комплексами ( + 1)-ого уровня (Рис.3.3) ( = 1, 2,...).

Формализуя эту идею, будем полагать, что комплекс двух молекул РНК иерархического уровня – это комплекс двух последовательностей с петлевыми участками, в которых энергия взаимодействующих мономеров перенормирована энергией петель иерархического уровня (+1). Пользуясь тем, что энергия петель в первом приближении пропорциональна длине (3.2), представим ее в виде:, где — длина петли, а — соответствующий -ому уровню коэффициент связывания. Подставляя в формулу (2.15) статистические веса петель,+ = /, получим выражение для определения свободной энергии комплекса двух случайных РНК–последовательностей 1 :

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

Будем считать, что -й и -й мономеры могут образовать связь, если:

а) участок [, 1] последовательности 1 не имеет связей с участком [, 1] подпоследовательности 2, где – минимальное количество нуклеотидов, необходимых для формирования петли определенного уровня (если < и/или <, то рассматриваются соответственно участки последовательностей [1, ] и/или по [1, ]);

б) 1-й мономер первой последовательности взаимодействует с 1-м мономером второй последовательности, и при замене ( 1), ( 1) выполняется а) (или б)).

В таблице 3.1 приведены значения для коэффициента связывания и минимальное количество нуклеотидов в петлях -ого уровня; вычисления проводились для случайных последовательностей равной длины = = 104. Длины последовательностей слабо Здесь, как и ранее, имеет смысл свободной энергии с обратным знаком влияют на средний коэффициент связывания, однако рассмотрение больших длин позволяет провести оценку для большего количества иерархических уровней. Отметим, что коэффициент связывания, определяемый по данной иерархической процедуре, медленно (логарифмически) стремится к 1 с ростом количества иерархических уровней (т.е.

при ). Логарифмическая зависимость обусловлена экспоненциальным ростом минимального числа мономеров, которые могут образовать петлю, = 3 + ( > 2) с увеличением номера иерархического уровня ( см. Табл. 3.1).

Таким образом, численно наблюдаемый коэффициент связывания (Рис. 3.2(а)) в действительности зависит от длин рассматриваемых последовательностей и полученное нами значение 0.92 лишь указывает на то, что последовательности длиной 400 1000 мономеров образуют структуру всего с двумя–тремя иерархическими уровнями.

Минимальная длина петли Коэффициент связывания Таблица 3.1 Вероятность связывания мономеров в зависимости от числа уровней в иерархической модели взаимодействия двух полимеров.

3.2 Распределение длин петель в РНК-подобных структурах Связывание двух РНК с петлевыми участками Было проанализировано распределение длин петель в структуре комплекса с петлевыми участками и внутрипетлевым взаимодействием. На Рис. 3.4 представлена зависимость () числа петель различной длины для структуры с петлевыми участками.

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

Действительно, величину = связывания мономера в структуре. Считая, что взаимодействие мономеров независимым, число петель длиной в структуре двух взаимодействующих сополимеров длиной можно оценить, как:

Такое распределение длин петель при 1 удовлетворяет очевидному соотношению =1 () = (1). Из Рис. 3.4 видно, что численные результаты хорошо аппроксимируются в логарифмическом масштабе прямой () =, где c хорошей точностью ln( 2 ) и ln(1 ) (см. (3.7)). Таким образом, в связывании сополимеров с петлевыми участками статистика петель выглядит в точности так, как происходит при независимом связывании мономеров. Однако стоит отметить, что модель независимого связывания дает хорошие результаты для последовательностей, в которых количество различных сортов мономеров 4. Для двухбуквенных и трехбуквенных алфавитов, взаимодействие сополимеров оказывается коррелированным, и формула (3.1) плохо описывает энергию оптимальной конфигурации.

Связывание двух РНК с внутрипетлевым взаимодействием Существенно иное поведение имеет статистика петель в комплексах с внутрипетлевым взаимодействием. На Рис.3.5(а) представлена зависимость числа петель с длиной по набору из 103 пар случайных последовательностей. Отметим особенности наблюдаемого распределения. Во-первых, для данной зависимости характерно степенное поведение. Показатель степенной зависимости для РНК разной длины меняется в интервале [1.38, 1.5]. Во-вторых, распределения для РНК с различной длиной совпадают, что позволяет проводить вычисления для набора коротких последовательностей.

В-третьих, при малых ( 5) характерно небольшое число петель с нечетной длиной и большое число петель с четной длиной. Последнее обстоятельство связано с тем, что для структуры комплекса с внутрипетлевым взаимодействием и = 0 характерно высокое значение средней энергии на один нуклеотид ( 0.92), обусловленное связыванием внутри петель, а образование петли малой длины с нечетным числом Рис. 3.4 Распределение длин петель в структуре комплекса с петлевыми участками.

Вычисления были выполнены для случайных последовательностей длины = 104, результаты усреднялись по набору из 105 сополимеров.

нуклеотидов приводит к потере, по крайней мере, одной возможной связи внутри петли. Таким образом, образование петель с нечетным числом мономеров энергетически невыгодно. Наконец, для распределения характерно наличие плато при больших, что обусловлено эффектом конечного размера (см., например, [86], где построена теория аналогичного эффекта).

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

Рис. 3.6). Каждому мономерному звену соответствует один шаг блуждания. Этот шаг направлен направо вверх, если мономер является «началом петли» (т.е. связан с мономером, расположенным после него по цепи), направо вниз, если он является «концом петли» (т.е. связан с мономером, расположенным до него по цепи) или горизонтально, если мономер не образует связи. Легко видеть, что такое построение задает соответствие между РНК-подобными вторичными структурами и так называемыми путями Моцкина [87] — состоящими из горизонтальных и диагональных участков дискретными случайными блужданиями в верхней полуплоскости, концы которых закреплены Рис. 3.5 (a) Распределение длин петель в структуре комплекса с внутрипетлевым взаимодействием. Вычисления проводились для последовательностей с длинами = = 75, 100 и 200, для каждого было выполнено 103 накоплений, для функция распределения сглаживалась по 10 соседним значениям); (б) Распределение путей Моцкина по длинам (длина пути случайного блуждания 200 шагов, количество накоплений -104, для 30 функция распределения сглаживалась по 10 соседним на оси абсцисс. Возвращение на ось абсцисс соответствует образованию одной петли в структуре комплекса. Как известно, [88], количество различных путей Моцкина (, ) длины с заданным количеством горизонтальных шагов определяется числами Каталана:

– биномиальные коэффициенты, ()/2 – числа Каталана. При (3.8) имеет асимптотическую зависимость (, ) 3/2 от длины пути. Было построено распределение длин петель для случайных путей Моцкина с вероятностью в численном моделировании значение вероятности образования связи, а вероятность горизонтального шага 1 2. Результат приведен на рисунке 3.5(б). Видно, что зависимость обладает всеми характерными свойствами, наблюдаемыми для распределения длин петель в структуре с внутрипетлевыми взаимодействиями.

Рис. 3.6 Вторичная структура РНК с пропусками и соответствующий ей путь Моцкина (а); Полностью связанная структура РНК без пропусков и соответствующий Представление структур РНК в виде путей Моцкина, статистика которых известна, позволяет сделать интересное наблюдение. А именно, показать, что для РНК-подобных структур характерно критическая зависимость структуры основного состояния в зависимости от числа различных сортов мономеров, используемых в последовательности.

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

4.1 Зависимость свободной энергии РНК-подобных структур от алфавита Рассмотрим случайную последовательность длиной и алфавита, образующую вторичную структуру типа РНК (Рис. 1.2(a)). Зададимся вопросом о том, к какому пределу стремится доля комплементарных пар в основном состоянии длинной ( ) цепи РНК. Другими словами, интерес представляет удельная (в расчете на одно звено) энергия основного состояния длинной РНК. Вначале, приведем доводы, подтверждающие наличие критического изменения удельной энергии в зависимости от алфавита.

Предположим, что существует критическое значение алфавита = такое, что при < доля связанных мономерных звеньев стремится к 1, тогда как при > предельная доля связанных звеньев меньше 1. Убедиться в этом можно следующим образом. Для того чтобы доля связанных звеньев в РНК-подобной структуре, образуемом случайной последовательностью, была равна 1, каждой последовательности из возможных должен соответствовать так называемый путь Дика (т.е. путь Моцкина, в котором нет горизонтальных шагов) (Рис. 3.6(б)). Количество путей Дика () длины определяется формулой (3.8) (() = (, 0)). При = 0 и при >> 1 () имеет асимптотическое выражение Заметим, что один и тот же путь Дика может описывать несколько РНК-подобных структур. Действительно, каждая пара подъем/спуск в пути Дика может быть, независимо от остальных, реализована разными способами (в случае РНК возможные варианты — это A-U, U-A, C-G и G-C). Таким образом, число различных первичных структур, для которых существуют полностью связанные вторичные структуры, не превышает Это оценка сверху, т.к., вообще говоря, одной и той же последовательности может соответствовать несколько различных РНК-подобных структур и, таким образом, несколько путей Дика. Тем не менее, естественно ожидать, что при число таких последовательностей с двумя и более полностью связанными вторичными структурами становится малым. В таком случае, сравнивая (4.2) с полным числом возможных первичных структур 0 (, ) =, можно записать ( 1):

Откуда, = 4. Подчеркнем, что несмотря на то, что эта оценка является грубой оценкой сверху, она демонстрирует характерное изменение свойств РНК-структур со случайной первичной структурой.

Таким образом, при < в пределе практически любой последовательности соответствует полностью связанная вторичная структура, и энергия оптимальной конфигурации на одну пару нуклеотидов стремится к 1, в то время как для случайных цепочек с > доля последовательностей, допускающих совершенную вторичную структуру, экспоненциально мала. Добавление горизонтальных шагов в пути случайных блужданий приводит к увеличению возможных РНК-подобных конфигураций (3.8), что позволяет сопоставить каждой случайной последовательности из ансамбля путь в случайном блуждании, соответствующий ее оптимальной вторичной структуре. Однако в этом случае доля связанных звеньев в оптимальной вторичной структуре остается в пределе меньше единицы. Путь Моцкина длиной, включающий горизонтальных шагов определяется (согласно (3.8)) как Для нечетных ( ) функция (, ) равна 0. Для четных ( ) воспользуемся (4.1) и формулой Стирлинга для оценки асимптотического поведения:

где введено обозначение = ( > 0). Последнее выражение (4.5) показывает рост для малых, но конечных.

Как много различных структур могут иметь один и тот же путь Моцкина (Рис. 3.6(а))?

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

А теперь оценим минимальное количество несвязанных мономеров (горизонтальных шагов в пути Моцкина), () = 1 (), в основном состоянии при > 4. Наибольшее количество структур, имеющих в основном состоянии долю несвязанных мономеров меньше или равной определяется выражением:

Для > 4 и = 0 эта сумма меньше 0 (, ) =, и растет с увеличением так, что при некотором величины (,, ) и 0 (, ) сравниваются. Для 1 сумму (4.7) можно оценить методом перевала. Введем обозначение где m = Для < m сумма в (4.8) определяется вкладом от верхней границы, тогда как для > m максимум достигается в точке m и, таким образом, не зависит от верхнего предела суммирования. Величина () определяется из уравнения (, ) = 0. На Рис. 4.1 представлена функция () = 1 (). Напомним, что данная оценка является верхней границей, так как не учитывает корреляции между оптимальными конфигурациями.

Оценка (4.8) сделана в предположении так называемого среднего поля: связывание на каждой паре подъем/спуск Рис. 3.6 происходит независимо с вероятностью 1/ и все пути случайных блужданий считаются статистически независимыми. В разделе 4.3 приводится более точная оценка критического алфавита, учитывающая корреляции между конфигурациями.

Результаты численного моделирования для РНК-подобных структур со случайной последовательностью звеньев различного алфавита представлены на Рис. 4.1. Для простоты предполагается, что комплементарные связи образуются согласно правилу A–A, т.е., только одинаковые мономеры могут комплементарно связываться. Напомним, что в реальных молекулах РНК действуют перекрестные правила комплементарности. Однако, анализ случайных последовательностей показал, что правила комплементарности незначительно влияют на свойства РНК-подобных структур. Тогда как, комплементарное связывания типа A–A позволяет исследовать цепочки не только с четным алфавитом как в случае перекрестного связывания, но и с нечетным. Соответственно, для каждого значения = 3, 4,..7 были построены зависимости удельной энергии =< > / от длины случайной первичной структуры. Как видно (Рис. 4.1(а)), удельная энергия при Рис. 4.1 (a) Зависимость удельной энергии от длины случайной последовательности с заданным алфавитом ; (б) зависимость предельного значения энергии от алфавита для последовательностей с дискретным алфавитом (красным), в модели Бернулли (синим), и верхняя оценка энергии (черным) в модели независимого связывания. Дополнительный график: зависимость предельного значения энергии от алфавита в модели Бернулли демонстрирует, что критический алфавит является ляется только функцией от (Рис. 4.1(б)). Результаты численного моделирования существенно расходятся с аналитической оценкой (Рис. 4.1(б)). Зависимость, полученная в численном моделировании имеет критическую точку топологического перехода = 2, которая является очевидной оценкой снизу. Действительно, рассмотрим произвольную двухбуквенную последовательность, например, и будем последовательно находить комплементарые пары (в предположении А–А связывания) следующим образом. Ближайшие соседи по цепи одного сорта образуют комплементарную пару, и далее, вычеркиваются из последовательности. Легко видеть, что такая процедура приводит к формированию РНК-подобной структуры. Рассматриваемая цепочка после первой итерации будет выглядеть:, последующее вычеркивание приведет к. Понятно, что данная процедура для любой случайной двухбуквенной последовательности приведет к тому, что, либо в остатке будет, либо последовательность будет полностью вычеркнута. Второй вариант означает, что все мономеры участвуют в формировании вторичной структуры, удельная энергия которой = 1. В случае остатка, данный участок цепочки образует конфигурацию с двумя пропусками, но, в термодинамическом пределе, = 1. Если структура образуется согласно перекрестным правилам комплементарности, то остатком будет |1 2 | букв одного сорта, где 1 и 2 — количество мономеров A и B соответственно. В случайной последовательности |1 2 | 1/. Таким образом, для случайных последовательностей с алфавитом = 2 можно записать:

Таким образом, аналитическое расcмотрение буквенных последовательностей позволило говорить о том, что критическое значение алфавита 2 4. Более строгое рассмотрение [89] показало, что критическое значение алфавита лежит в интервале:

Подводя итог, подчеркнем еще раз, что при изменении алфавита, используемого в первичной структуре случайной последовательности, существует переход от полностью связанной РНК-подобной структуры до структуры с конечной долей несвязанных мономеров. Такой переход в работе называется топологическим. Критическая точка топологического перехода принадлежит интервалу (2 < < 3), т.е. эффективно является нецелым.

Как можно трактовать нецелый алфавит в случайных последовательностях типа РНК? Далее, попробуем ответить на этот вопрос.

4.2 Топологический переход в модели Бернулли Модель случайной последовательности с эффективно нецелым алфавитом может быть построена следующим образом. Будем считать, что матрица контактов в уравнении (2.20) является случайной: вероятность того, что, = 1, равна, а вероятность, = 0 равна соответственно 1. То есть теперь случайная последовательность характеризуется не первичной структурой — последовательностью мономеров из различных типов, как это было раньше, а некой матрицей контактов, (, )-элемент которой Для доказательства использовалось понятие максимального паросочетания без пересечений на случайном слове, работа была выполнена после доклада в ИППИ РАН (май 2012) разрешает или запрещает образование комплементарной пары между и мономером цепи. Мономеры цепи в данной модели не различаются по сортам и, в целом, любой мономер может образовать связь с любым другим в цепи, однако, в среднем, вероятность такого события равна. Каждой последовательности в рассматриваемой модели можно сопоставить граф Эрдёша–Реньи, изображающего все возможные контакты между мономерами. Основное отличие данной модели от дискретных буквенных последовательностей — нарушение свойства транзитивности. Если 1-й мономер может образовать связь со 2-м, а 2-й с 3-м, отсюда, вообще говоря, не следует (как это было для последовательностей с дискретным алфавитом), что 1-й мономер может связаться с 3-м. Однако, как, например, уже упоминалось, подобная модель бернуллиевского сравнения в задачах выравнивания случайных последовательностей является хорошей аппроксимацией.

Вероятности случайной матрицы контактов соответствует алфавит, равный:

Таким образом, оказывается возможным генерировать случайную последовательность с любым нецелым значением алфавита. На Рис. 4.1(б) приведена зависимость удельной энергии в термодинамическом пределе от алфавита (4.10), полученная в численном моделировании. Во-первых, отметим, что значения для бернуллиевского алфавита не более, чем на 1% отличается от соответствующих величин для случайных последовательностей с дискретным алфавитом, что оправдывает применимость данной модели. Случайный бернуллиевский полимер характеризуется критической вероятностью. Для >, в термодинамическом пределе, = 1 (так называемая «полочка»

на зависимости удельной энергии (см. дополнительный график на Рис. 4.1(б)), что соответствует полностью связанной вторичной структуре, тогда как для <, даже в пределе бесконечной длины, основное состояние характеризуется () количеством несвязанных мономеров. Критическое значение вероятности согласно (4.10) соответствует критическому значению алфавита = 2.6. Таким образом, модель Бернулли позволяет численно получить точку перехода.

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

Рис. 4.2). Естественно ожидать, что в пределе (Рис. 4.2(а)), функция () вырождается в ступенчатую функцию. Скейлинг-анализ полученных зависимостей () обеспечивает критическое значение = 0.37, что соответствует алфавиту:

Рис. 4.2 Зависимость доли полностью связанных РНК-подобных структур в ансамбле случайных первичных структур различной длины (а) от параметра модели Бернулли;

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

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

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

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

где 1 и 2 — некоторые константы. Для допереходной фазы ( > ) характерно экспоненциальное приближение к предельному значению удельной энергии ( = 1), тогда как в области больших алфавитов ( < ) энергия приближается к своему предельному значению степенным образом (Рис. 4.3). Показатель степени в (4.11) находится в пределах [0.75, 1] (сравните с (Рис. 4.1(а))). В допереходной области случайная последовательность из алфавита может быть охарактеризована некоторой релаксационной длиной, указывающей на характерный масштаб длин, на котором энергия основного состояния сходится к своему предельному значению = 1. Ясно, что зависимость релаксационной длины от вероятности имеет вертикальную асимптоту в точке =.

Естественно ожидать, что асимптотическое поведение () зависит от выбранной модели случайного полимера, в частности от правил комплементарности — см. (4.9).

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

Рис. 4.3 Асимптотическое поведение удельной энергии () до (а) и после (б) топологического перехода. Зависимость ( ()) в логарифмическом масштабе (а) и двойном логарифмическом масштабе (б) (см. (4.11)).

4.3 Аналитическая оценка критической точки топологического перехода в модели Бернулли 4.3.1 Метод среднего поля Для простоты переформулируем задачу в терминах планарных диаграмм (Рис. 1.7).

Рассмотрим граф, вершины которого (мономеры вдоль цепочки) перенумерованы, а матрица контактов — матрица инцидентности графа. Задача о полностью связанной РНК-подобной структуре на данном графе сводится к вопросу о том, как выбрать среди разрешенных контактов /2 связей, которые обеспечивают планарную структуру на заданном случайном графе, т.е. все вершины входят в конфигурацию ровно один раз и любые пары связей (1, 1 и (2, 2 ) удовлетворяют соотношению [91]:

Другими словами, как разместить /2 непересекающихся арок, принимая во внимание ограничения, накладываемые матрицей. В модели Бернулли каждый элемент равен 1 либо 0 с соответствующими вероятностями и 1, кроме того, матрица контактов — симметричная с нулевыми диагональными элементами:

Здесь () и () — дельта-функция Дирака и функция Хевисайда, соответственно. Для = 1 (когда все элементы равны 1), количество всех возможных арочных структур, удовлетворяющих (4.12) определяется числами Каталана (см. (4.4)) Когда = 1, некоторые из конфигураций # запрещены матрицей контактов. Введем обозначение 1 — вероятность, того, что одна выбранная из # конфигурация разрешена.

Очевидно, что Аналогично, определим как вероятность, что диаграмм из # разрешены, для = 2, например где 12 2 равно количеству общих арок для двух случайно выбранных планарных диаграмм, усредненному по ансамблю #. Для 3 можно записать:

Величины могут быть вычислены с любой точностью. К примеру, 2 лежит строго в интервале [1/15, 1/14.8]. Вероятность иметь по крайней мере одну планарную конфигурацию для данной заполненности матрицы (4.13) определяется как:

Предполагая, что все диаграммы в ансамбле # независимы, т.е., =, для из (4.18) можно записать:

В пределе больших, величина равна либо нулю, либо единице, в зависимости от соотношения между # и 1. Используя (4.15), для критического значения вероятности можно записать уравнение:

Условие (4.20) можно интерпретировать как то, что переход наблюдается в точке, при которой плотность единиц в матрице контактов такая, что в среднем разрешена только одна планарная конфигурация. Вспоминая, асимптотику чисел Каталана (4.16), для критического значения вероятности получим = 1/4, что совпадает с верхней оценкой = 4 из (4.2).

4.3.2 Комбинаторная оценка Предположение о независимости планарных конфигураций соответствует так называемому приближению среднего поля. Естественным следующим шагом является введение ненулевых корреляций между конфигурациями: = 0. Чтобы учесть корреляции между различными планарными диаграммами, поступим следующим образом.

Перепишем (4.20) как:

где () — некоторая функция, учитывающая корреляции между планарными диаграммами. Основная идея дальнейшего рассмотрения следующая: арки разной длины встречаются в оптимальной планарной конфигурации с различной вероятностью. Рассмотрим полностью связанную планарную конфигурацию, состоящую из = арок, соединяющих точек. Возвращаясь к представлению планарных диаграмм через пути Дика (см.Рис. 3.6), можно увидеть, что арка между -ой и -ой точками возможна, только если -й и -й шаг имеют одну и ту же пространственную координату. Тогда можно определить вероятность арки между -ой и -ой точками как:

В знаменателе правой части (4.22) стоит суммарное число возможных шагов вверх/вниз на длине ( + 1), в числителе — “1” соответствуют выбору шага вверх и вниз на позициях и соответственно; число Каталана (1)/2 описывает все возможные конфигурации петли между парой (, ) (так как -й и -й шаги находятся на одной высоте, петля между ними должна быть тоже путем Дика). Вероятности (, ) зависят только от длины арки ( ) и не равны нулю только для арок нечетной длины, т.е., вероятность того, что в -ой позиции находится левая граница арки (шаг вверх).

Отметим, что доля коротких арок чрезвычайно высока. Действительно, вероятность, в типичной арочной конфигурации иметь арку длиной = 1 равна 1, арку длиной = 3, уже, и т.д.. С другой стороны, количество всех возможных кратчайших арок — ( 1). Поэтому, в типичной конфигурации среди них должны быть «разрешены».

Естественно, что веса таких коротких арок в бернуллиевской модели (элементы,+ матрицы контактов) выше, чем длинных арок.

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

Так как общее число длинных арок порядка 2, будем считать, что длинные арки выбираются независимо друг от друга с вероятностью. И, таким образом, вклад от длинных арок в функцию () равен /4.

Иная ситуация при выборе кратчайших арок длиной “1”. Для бернуллиевского полимера с матрицей контактов только единичных арок разрешены. Таким образом, выбор коротких арок для оптимальной конфигурации без пропусков оказывается сильно ограниченным. Вероятность выбрать /4 непересекающихся арок из разрешенных непересекающихся единичных арок из всех ( 1) возможных (Рис. 4.4). Единичные арки можно рассматривать как стенки ящиков, тогда задачу можно переформулировать следующим образом. Будем интересоваться количеством способов, которыми можно заполнить ( 1) ящика /2 свободными точками (шарами). Результат известен из комбинаторики и = / Можно считать, что среди них 3 1 арок разрешены первичной структурой полностью связанной РНК-подобной структуре случайного полимера. Учет корреляций между планарными конфигурациями на уровне единичных дуг приводит к следующими выражению для () (4.21):

В пределе, после упрощений, получим:

случайного графа ( свободных вершин) аналогичен комбинаторной задаче о Подставляя этот результат в (4.21):

4.3.3 Матричный подход Еще один подход оценки критического алфавита основан на матричном описании вторичной структуры РНК (1.5). Напомним, что статистическую сумму (, ) случайного полимера можно представить через случайные эрмитовы матрицы как (см. (1.7),(1.8)):

где В отсутствии замороженного беспорядка, т.е., если 1, задача (4.26) может быть решена точно. В частности, множитель /2,0 перед /2, описывающий планарные конфигурации с /2 арками (1.10), т.е. полностью связанные структуры, вносит наибольший вклад в общую статистическую сумму полимера и определяется числами Каталана:

Как и ранее будем вычислять функцию () в (4.21) усредняя статистическую сумму (, ) по распределению (4.13). Для этого выполним стандартное преобразование Хаббарда-Стратоновича и будем интегрировать по с весом (4.13):

Величина 0 соответствует единичной матрице контактов с дополнительным фактором. Учет только этого слагаемого, после обратного преобразования ХаббардаСтратоновича, приводит к () =, и оценке = 4, совпадающей с оценкой в предположении среднего поля. Действие (4.31) сдвигает значение в сторону меньших значений. Но, так как содержит бесконечное число слагаемых (4.31), теория возмущений в данном случае неприменима. В этой связи было предложено следующее приближение: все поля { }=1,..., в (4.31) эквивалентны, поэтому можно считать, что в среднем, tr( )0 не зависит от (, ). В рамках данного средне-полевого приближения можно сделать замену = 0, где:

Упрощение выражения (4.32) приводит к следующему уравнению на пропагатор :

Выражение (4.33) дает = 2 log 1 11/, и окончательно можно написать:

Подстановка (4.35) в (4.21) приводит к оценке критического алфавита * = 0.4551.

Большая расходимость полученного результата с численным = 0.37 означает, что предложенного приближения недостаточно для описания топологического перехода.

4.4 Переход случайной РНК в замороженное состояние, ограниченный топологическим переходом Рассмотрим как данный топологический переход ограничивает фазовый переход в замороженное состояние 1.4. Отметим, что аналогичный вопрос исследуется и в теории перколяции, где тоже предполагается взаимосвязь перколяционного перехода и температурного фазового перехода, наблюдаемого, например, в модели Изинга [92].

Были проанализированы температурные зависимости свободной энергии пинча (2.2) случайной последовательности в модели Бернулли разной вероятности. Как уже обсуждалось, температура перехода в замороженное состояние непосредственно связано со средним числом пропусков в структуре основного состояния. В [38] было показано, что температура перехода не превосходит * где — среднее число пропусков на пару мономеров, а определяется из зависимости наибольшего общего непрерывного сегмента двух половинок последовательности РНК: = 1 ln (см. Рис. ??). Известно, что для цепочек РНК = ln 2. Для случайного бернуллиевского процесса определяется как = ln(1/) [79]. Таким образом, выражение (4.36) можно переписать в виде Доля несвязанных мономеров растет с ростом алфавита 1/ сильнее, чем логарифм (см. Рис. 4.1(б)) и из (4.37) непосредственно следует, что в допереходной области ( > ) фазовый переход в замороженное состояние наблюдаться не будет. Температура перехода эффективно равна нулю, т.е., случайный полимер во всем температурном диапазоне находится в расплавленной фазе. Данное предположение дополнительно подтверждается наблюдением того, что для случайных последовательностей с алфавитом = 2 переход имеет место только при накладывании ограничений на структуру, а именно, введением минимального размера петли [39].

Результаты численного моделирования представлены на Рис. 4.5. Был проанализирован температурный коэффициент ( ) (2.4) для последовательностей с разной вероятностью. Температура перехода определяется точкой, в которой нарушается линейная зависимость ( ) = 2, характерная для расплавленной фазы. Из полученных данных видно, что температура перехода уменьшается с ростом вероятности и в допереходной области становится равной нулю ( = 0.5 на Рис. 4.5). Вблизи критического значения численный эксперимент усложняется тем, что корректный анализ требует рассмотрения достаточно длинных случайных цепочек (с длиной, превышающей соответствующую релаксационную длину (), см. (4.11)), что приводит к существенному увеличению времени численного моделирования. Также стоит отметить, что в связи с наблюдаемой степенной зависимостью свободной энергии основного состояния от длины последовательности ((4.11)), аппроксимация уравнением (2.2) вблизи точки = 0, вообще говоря, неверна.

Предполагается, что критическая точка топологического перехода между полностью связанной РНК-подобной структурой и структурой с пропусками является пороговым значением для термодинамического перехода. В области последовательностей > возможна только расплавленная фаза вне зависимости от температуры. Рис. 4.6 покаa(T) Рис. 4.5 Зависимость коэффициента ( ) (1.4) для случайной последовательности зывает фазовую диаграмму на (, ) плоскости. Это предположение подтверждается исследованием энергии пинча от длины случайной последовательности в точке = 0.

Точка пересечения зависимостей для разных длин (см. дополнительный график на Рис. 4.6) разделяет два топологических режима и близка к наблюдаемому критическому алфавиту.

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