Monthly Archives: April 2012

Секвенирование методом RAD. Подготовка библиотеки

Секвенирование методом RAD (Restriction site Associated DNA) – один из наиболее эффективных методов секвенирования следующего поколония (next-generation sequencing), особенно выгодный в условиях, когда есть необходимость обработать данные для большого числа организмов (например, при изучении эволюционной динамики популяций) или собрать новый геном (de-novo genome assembly) без референсного генома. Для того, чтобы дать возможность читателям ознакомиться подробнее с этой технологией, мы разбиваем  тему на отдельные блоки:

  • Блок 1. Подготовка библиотек и секвенирование
  • Блок 2. Алгоритмы обработки сырых данных (reads)
  • Блок 3. Практическое использование обработанных данных для научных задач

Блок 1. Подготовка библиотек и секвенирование

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

В основе метода RAD лежит простая идея (описывется для секвенирования с помощью платформы Illumina, схематическое изображение процесса представлено на картинке ниже):

Нуклеотидные последовательности в ДНК-образце для определенного организма  разбиваются на блоки с помощью рестриктразы (restriction enzyme  – например, EcoRI или SbfI). Затем присоеденяется первый адаптор Р1 (P1 adaptor – например, Solexa от Illumina или его варианты), содержащий уникальную последовательность ДНК – т.н. баркод (barcode). Этот баркод позволяет уникальным образом идентифицировать каждый ДНК-образец. Образцы смещиваются (pooled, “уникальность” потом можно будет определить по баркодам – отличающимся для каждого организма), нуклеотидные последовательности разбиваются еще раз на более мелкие кусочки с помощью механического воздействия (shearing) и с помощью гелевого электрофореза определяются блоки нуклеотидов определенной длины и изолируются (для платформы Illumina длина может быть 300-700 bp). Обратите внимание, что после механического воздействия часть нуклеотидных блоков будет содержать первый адаптор, а часть – нет. Затем присоединяется второй адаптор P2 c так называемыми расходящимися концами (diverging ends). Затем, с помощью праймеров, уникальных к адапторам P1 и P2 (primers P1-/P2-specific), амплифицируются только фрагменты, содержащие оба адаптора. Это достигается за счет того, что адаптор Р2 имеет расходящиеся концы, а такой конец не может быть присоединен к праймеру Р2, до тех пор, пока не будет “закончен” при амплификации первого адаптора Р1. Таким образом на выходе мы получаем следующую структуру:

адаптор Р1 – баркод – участок, к которому присоединяется рестриктраза (например, для SbfI TGCAGG) – ХХХ-ХХХХ bp последовательность ДНК – адаптор 2

Из полученной смеси опять выделяют фрагменты определенной длины (примерно 200-500 bp) и именно это представляет из себя законченную библиотеку которую потом можно секвенировать на платформе Illumina. “Считывание” на платформе Иллюмина будет происходить начиная с баркода, обеспечивая таким образом возможность уникальной идентификации последовательности ДНК по отношению к исходным образцам. Детальный протокол (для подготовки библиотеки в собственной лаборатории) можно прочитать вот здесь.

Из очевидного: не забудьте, что вам нужно выделить ДНК перед тем, как проводить подготовку библиотеки.

 

Rad sequencing - library preparation

Rapid SNP Discovery and Genetic Mapping Using Sequenced RAD Markers: http://www.plosone.org/article/info:doi/10.1371/journal.pone.0003376

Некоторые советы полученные от разработчиков метода и других исследователей:

Подготовка баркодов.

  1. Уникальные последовательности должны уникально различимы (!), то есть с учетом того, что могут происходить ошибки считывания у платформы Illumina. Например, если вы выбираете  баркод вида ATTGT для образца 1 и баркод вида ATTCT для образца 2, есть вероятность получения образца с ошибкой чтения именно в области замены G’C – таким образом, отличить образцы вы не сможете и будете считать образец 1 образцом 2. Поэтому рекомендуется брать разницу хотя бы в 2, а лучше в 3 нуклеотида между баркодами.
  2. Если есть желание делать баркоды разной длины – это тоже может быть эффективным иструментом для разделения образцов и повышения качества reads recovery.
  3. Не используйте повторяющиеся нуклеотиды для создания баркодов (например AAAAA или AAAGG). Последовательность должна быть достаточно вариабельная.

Покрытие.

Очень важный аспект – расчет конечного покрытия отсеквенированных последовательностей. То есть сколько раз будет отсеквенирован блок выше. Почему это важно? потому что после проведения анализа будет необходимо отделить истинные нуклеотидные замены от ошибок секвенирования. Для того, чтобы удалить ошибки (которые будут распределены случайно), нужно иметь повторные реплики каждого секвенируемого блока – то есть в англоязычной терминологии покрытие блоков. Например, 25, 30, 45 раз. Тогда статистически будет проще узнать какие замены являются истинными, а какие – нет. Расчет покрытия зависит от:

  • числа образцов
  • длины рестриктразы
  • размера генома

Найти расчетные формулы для примерной оценки и обсуждение можно вот здесь. Рекомендуемое покрытие – 25-50 раз.

Стоимость.

Стоимость анализа будет складываться из следующих составляющих:

  • получение образов и выделение ДНК (в пределах 1000 долларов, если есть своя лаборатория – дешевле)
  • синтез адапторов с баркодами (расчет на 2009 год) – на 50 адапторов примерно 8 000 долларов при объеме 0,2 umol
  • секвенирование на машине Illumina (зависит от лаборатории и страны – но цены можно найти в интернете)

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

Скрипт для быстрого анализа покрытия генов в базе GENBANKBiopython script to analyse gene coverage in GENBANK for a given taxa

В связи с тем, что объемы информации в базе генетических данных GENBANK растут невероятными темпами, очень часто пользователям приходиться выяснять покрытие генов для определенной таксономической группы. Например, перед тем как провести реконструкцию филогенетического дерева, хочется узнать какие гены являются наиболее полно представленными в базе. При “ручном” решении задача является очень трудоемкой – так как нужно отсмотреть сотни или даже тысячи страниц с результатами поиска вида через интерфейс NCBI, выбрать наиболее представленные гены, определиться с перспективными референсными последовательностями и лишь потом – загрузить сами данные из базы.

Для того, чтобы упростить эту задачу был написан скрипт на biopython, который позволяет автоматически собрать статистику по представленности генов в базе данных GENBANK  в пределах конкретной таксономической группы. В качестве единственного параметра скрипт принимает название таксономической группы (например, вид, род, семейство, порядок, и т.д.). Важно понимать, что скрипт будет загружать на локальный компьютер данные в формате genbank перед тем, как проводить парсинг записей и сбор статистики, поэтому не рекомендуется запращивать, к примеру, все набор данных для Angiosperms – загруженный файл будет огромным по размерам, а ваш IP может быть заблокирован NCBI.

Картинка для featured image поста взята отсюда

Cкрипт ниже написанный на R автоматически обрабатывает выходной файл скрипта на питоне и позволяет сразу получить таблицу, в которой первая колонка будет уникальное название гена, затем сколько раз ген был встречен в данной таксономической группе, ID референсной генетической последовательности (наиболее длинной последовательности) и число уникальных видов для которых был отскевенирован данный ген.

 В связи с тем, что объемы информации в базе генетических данных GENBANK растут невероятными темпами, очень часто пользователям приходиться выяснять покрытие генов для определенной таксономической группы. Например, перед тем как провести реконструкцию филогенетического дерева, хочется узнать какие гены являются наиболее полно представленными в базе. При “ручном” решении задача является очень трудоемкой – так как нужно отсмотреть сотни или даже тысячи страниц с результатами поиска вида через интерфейс NCBI, выбрать наиболее представленные гены, определиться с перспективными референсными последовательностями и лишь потом – загрузить сами данные из базы.

Для того, чтобы упростить эту задачу был написан скрипт на biopython, который позволяет автоматически собрать статистику по представленности генов в базе данных GENBANK  в пределах конкретной таксономической группы. В качестве единственного параметра скрипт принимает название таксономической группы (например, вид, род, семейство, порядок, и т.д.). Важно понимать, что скрипт будет загружать на локальный компьютер данные в формате genbank перед тем, как проводить парсинг записей и сбор статистики, поэтому не рекомендуется запращивать, к примеру, все набор данных для Angiosperms – загруженный файл будет огромным по размерам, а ваш IP может быть заблокирован NCBI.

Картинка для featured image поста взята отсюда

 

Cкрипт ниже написанный на R автоматически обрабатывает выходной файл скрипта на питоне и позволяет сразу получить таблицу, в которой первая колонка будет уникальное название гена, затем сколько раз ген был встречен в данной таксономической группе, ID референсной генетической последовательности (наиболее длинной последовательности) и число уникальных видов для которых был отскевенирован данный ген.