Monthly Archives: July 2015

Загрузка данных из NCBI для реконструкции больших филогений – часть 2

Публикуем продолжение (вторую часть) туториала по получению данных для реконструкции больших филогенетических деревьев. Этот туториал пока доступен на английском языке на сайте InsideDNA.me, но скоро появится и тут на русском. Второй инструмент для загрузки больших объемов данных из NCBI – geneCoverage2fasta позвляет быстро получить фаста-файлы для всех референсных генов  с помощью BLAST. В английском туториале приводим объяснение почему нужно использовать BLAST, а не, к примеру, названия генов или видов для парсинга базы NCBI.

Анализ ортологов в геномах трех бактерий Ralstonia с помощью OrthoMCL

Один из наиболее мощных инструментов для кластеризации ортологов – это программа OrthoMCL. Несмотря на свою популярность, это также одна из наиболее трудно устанавливаемых и настраиваемых программ в биоинформатике. Мы постарались максимально упростить работу с этим инструментом в рамках InsideDNA и описываем ниже полный пайплайн для анализа ортологов в OrthoMCL / MCL .

Использование сайта InsideDNA бесплатно и открыто для всех!

1. Входные данные

Мы будем обрабатывать геномы для трех видов бактерий рода Ralstonia полученных с ftp NCBI. Все три генома можно скачать в одном общем архиве вот здесь (792Kb). Исходные данные для OrthoMCL могу быть представлены как в виде fasta файлов с аминокислотными последовательностями, так и нуклеотидными. Кроме того, это могут быть как аннотированные данные, так и просто результаты секвенирования RNA-seq или single reads/paired end данных Illumina. Чаще всего анализируются данные по нескольким геномам.

2. Сборка всех нужных инструментов в проект

Залогиньтесь в приложение InsideDNA и перейдите во вкладку Tools. На первом этапе мы создадим проект Bacterial_orthomcl

 1

и добавим туда все необходимые для анализа инструменты (7):  mcl,  orthomclAdjustFasta, orthomclFilterFasta, orthomclBlastParser, orthomclLoadBlast, blastall,  formatdb.

2

3. Запуск orthomclAdjustFasta

Первый этап  – это подготовка входных fasta-файлов. На этом этапе мы форматируем строку названия последовательности. Например, в нашем случае мы приводим заголовки каждого из трех входных fasta файлов из вида:

к виду:

>RE|NP_942643.1

Где RE – это краткое название вида (например, Ralstonia_eutropha) Процедуру коррекции заголовков нужно повторить для всех входных fasta файлов (например, для всех видов бактерий, в которых мы хотим найти ортологи). Мы используем вот такие параметры для вида Ralstonia_solanacearum:

 3

В нашем случае, поле с уникальным номер последовательности находится на четвертой позиции – но в случае других данных, позиция может отличаться. Нам нужно повторить процедуру для всех трех геномов

 4

В результате мы получим три откорректированных файла

 5

  Вот так выглядит fasta-файл до коррекции

 6

И после

 7

Переместим данные в отдельную папку adjusted

 8

4. Запуск orthomclFilterFasta

На следующем шаге отфильтруем все последовательности во всех геномах по длине и числу стоп-кодонов (минимальная длина протеина – 20 АА и максимальный процент стоп-кодонов – 20%) и объеденим их в один общий файл. 9

Именно для этой команды мы перемещали все откорректированные геномы в отдельную папку.

10

На выходе мы получим единый файл со всеми данными из всех трех геномов который называется goodProteins.fasta. Также мы получим файл poorProteins.fasta с отфильтрованными (плохими) последовательностями.

 11

5. Запуск formatdb

Создадим новую папку db

 12

И переместим файлы goodProteins.fasta и poorProteins.fasta в эту папку

 13

Теперь запустим форматирование данных для reciprocal blast с помощью команды formatdb. Эта команда трансформирует фаста-файл в специальный формат базы данных BLAST.

 14

На выходе мы получаем три дополнительных файла – phr, pin, psd

 15

5. Запуск Blastall

Теперь запустим одну из наиболее важных команд, которая и выполняет reciprocal blast (one-to-one blast – бласт один к одиному) и проверяет парное сходство между всеми последовательностями во всех геномах. Откроем окно параметров Blastall. Во-первых, нам нужно указать нашу пользовательскую базу BLAST – ту самую, что мы форматировали с formatdb. Приставка этой базы – goodProteins.fasta Во-вторых, мы должны передать серверу все входные файлы относящиеся к базе – phr, pin, psd, fasta – выбираем их в третьем поле сверху В-третьих, мы задаем т.н. query file – это собственно говоря тот же самый goodProteins.fasta – именно потому, что мы шлем запрос с данными к базе данных на основе этого же файла процедура называется reciprocal blast.

 16

В-четвертых, мы задаем выходной файл – blastallvall_prot.csv В-пятых, наши данные – это аминокислотные последовательности – следовательно мы используем алгоритм Blastp В-шестых, мы можем указать e-value – пороговое значение по которому будут определяться парные последовательности. В-седьмых, важно не забыть указать формат выходного файла – нам нужен тип 8 – tabulated format. Именно он может быть корректно прочитан программой orthoMCLBlastParser на следующем этапе.

 17

Для blastall лучше указать очередь помощнее (8 ядер и 30 Gb RAM)

 18

В результате обработки программы будет создан blastallvall_prot.csv с содержимым вида:

 19

6. Запуск orthoMCLBlastParser

После того, как мы получили файл blastallvall_prot.csv мы должны его обработать программой OrthoMCLBlastParser. Для этой программы нам потребуется как файл blastallvall_prot.csv, так и директория adjusted с тремя откорректированными фаста-файлами с геномами бактерий. Назовем выходной файл similarSequences.txt.

220

Это программа в результате создаст файл similarSequences.txt в папке blast_parser со следующим содержимым

221

7. Запуск orthoMCLLoadBlast

Теперь мы готовы к финальному шагу с использованием OrthoMCLLoadBlast – на вход нам надо предоставить файл similarSequences.txt. На вход этой программе нам надо передать всего лишь один файл, а также указать название лог-файла. Остальные выходные файлы будут автоматически созданы в той же директории, где расположен файл similarSequences.txt (то есть в нашем случае –  в папке blast_parser)

222

В результате работы этой программы будут созданы следующие выходные файлы и папки

223

Внутри папки pairs находятся следующие файлы

224

Нас же интересует текстовый файл mclInput – который и требуется для финально кластеризации ортологов с помощью программы MCL

225

8. Запуск MCL

Финальным шагом анализа ортологов в OrthoMCL является запуск программы MCL. Для кластеризации нам необходимо указать созданный mclInput , а также выбрать пороговое значение для кластеризации на графе – параметр -I, который мы ставим равным 1.5. Важно понимать, что с подбором данного параметра нужно экспериментировать и его оптимальное значение будет определяться тем, насколько эволюционно дивергентными являеются виды.

226

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

227

Всего, в случае с тремя геномами бактерий Ralstonia мы получили 562 группы ортологов. Обычно т.н. conserved genes – то есть гены представляющие core genome – представлены однократно в каждом геноме и расположены в нижней части csv файла, в то время как псевдогены (к примеру, ретротранспозоны) будут расположеные в верхних строчках csv файла и будут представлены общирными группами с сотнями генов. И те и другие группы ортологов интересны – например, core genome genes могут быть использованы для реконсрукции консенсусного дерева, в то время как ретротранспозоны – для анализа истории дупликаций или горизонтального переноса генов.

Оценка покрытия генов и секвенированных видов в NCBI для реконструкции больших филогений

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

Для простоты рассмотрим следующий характерный случай – вы отсеквенировали несколько видов из семейства Balanophoraceae для двух участков: 5.8S ribosomal RNA и 18S ribosomal RNA и теперь хотите:

  • Оценить какие виды Balanophoraceae были отсеквенированы ранее и присутствуют в базе Генбанк
  • Загрузить последовательности для данных генов
  • Объединить последовательности с тем, что вы получили в лаборатории
  • Реконструировать филогенетическое дерево из всех доступных генов для всего семейства Balanophoraceae

Чтобы не делать этот туториал слишком длиным, мы рассмотрим первый шаг здесь, а 2, 3 и 4 в следующих туториалах (подпишитесь на рассылку от InsideDNA.me)

Оценка того, какие виды Balanophoraceae species были отсеквенированы и для каких генов

  1. Регистрация.

Создайте аккаунт на сайте InsideDNA. После того, как вам придет письмо-уведомление, нужно будет заполнить небольшую анкету (~3 минуты) и после этого вы получите полный доступ на сайт.

1

  1. Общий принцип работы InsideDNA

В InsideDNA есть три основые закладки для работы – они расположены в главном меню навигации и называются Tools, Tasks и Files. В этих закладках можно делать следующее:

  • В меню Инструменты (Tools) можно искать различные биоинформатические и биологические инструменты (скрипты, программы, пайплайны), организовывать их в проекты и запускать
  • В меню Задачи (Tasks) можно мониторить процесс выполнения задачек и совершенать над ними некоторые операции (ниже об этом подробнее)
  • В меню Файлы (Files) можно управлять своими данными и файлами так, как если бы они были расположены на вашем компьютере. На данном этапе каждый пользователь получает 10 Гигобайт места сайте.
  • 2

     

    1. Создание нового проекта.

    Для того, чтобы создать новый рабочий проект нажмите кнопку +Add new project в левом меню 3   Назовите проект Balanophoraceae.

    4

    Теперь, мы можем добавить инструменты в проект – это удобно для работы, так как позволяет легко организовывать инструменты и не тратить время на их поиск и в будущем настройку. Добавим два инструмента в проект – в строке поиска наберите: geneCoverage

    5

    Затем нажмите на кнопку плюс на инструменте и выберите из списка проект Balanophoraceae.

    6

    Повторите операцию для инструмента geneCoverage2fasta – добавьте его в проект.

    1. Инициализация и настройка задачи

    Далее, просто выберите проект Balanophoraceae в левом меню – там вы увидите два добавленных инструмента и нажмите на кнопку Run Tool на geneCoverage. Более детальную информацию о каждом инструменте можно получить нажав кнопку Read More.

    7

     

    1. Выберите настройки инструмента

    После того, как вы нажмете на Run Tool появится окно настройки инструмента. 8 Здесь нужно указать: название задачи, и параметры инструмента. Затем мы нажмем окно препросмотра задачи и запустим ее. Укажите такое имя задачи, которую вам легко запомнить и идентифицировать в будущем – например, Balan_geneCov 9 geneConverage очень простой инструмент. Вам нужно всего лишь указать таксон, в пределеах которого будет осуществляться оценка покрытия по базе ГенБанк и выбрать директорию в которую будут сохранены результаты. Мы указываем следующее:

    10

    Balanophoraceae в качестве входного таксона (т.е. мы ищем все последовательности для видов внутри данного семейства).
    Далее нам нужно нажать кнопку Browse, которая вызывает диалоговое окно мини Файл менеджера. Этот мини файл менеджер позволяет быстро указывать входные и выходные файлы и папки для инструмента, а также создавать новые папки. Так как у нас пока не создано папки проекта нажмен на кнопку Плюс и создадим новую директорию Balanophoraceae_project. Это будет наша рабочая директория для данных

    11

    Выберем эту папку (нажмите Select). Все выходные данные автоматически будут скопированы в эту папку.

    13

    Хорошая практика – всегда проверять синтаксис команды перед отправкой – для большинства инструментов это позволяет еще и познакомиться с тем, какие параметры мы бы могли указать в командной строке. Для предпросмотра задачки нажмите Preview task

    14

    Если вы довольны выбраными настройками, последний важный шаг это выбор Очереди для задачи. Разные Очереди дают доступ к разным вычислительным мощностям и хорошая идея сначала запустить задачу в более «слабой» очереди преждем чем переходить к мощным серверам. Для нашей задачи нам хватит очереди с 4Гб рам.

    15

      После того, как вы выбрали очередь нажмите кнопку Submit task. После отправки перед вами появится диалоговое окно, где вам будет предложено либо остаться на текущей странице с Настройками инструмента (Stay on the page), либо перейти в раздер мониторинга задач (View submitted task). Остаться на текущей странице бывает удобно, когда необходимо отправить много похожих задачек меняя например только входные данные или некоторые параметры в настройках. Выберите Go to Tasks

    16

    1. Мониторинг задачи.

    В Задачах есть несколько опций: 1) Задачи сгруппированы по статусам. Обычно, если задача обрабатывается без ошибок, то тогда в начале она появится в группе со статусом Running и после завершения будет перемещена в группу Completed.

    17

    Иногда, если вы отправили более 5 задач одновременно, то 6ая задача и далее будет помещена в группу Suspended (ожидание). Как только освободиться слот в текущем блоке задачек, одна из задач в группе Suspended будет перемещена в группу Running. В некоторых случаях, задачи могут ломаться и тогда они появляются в группе Failed. Основные причины по которым задачи могут ломаться – это недостаток вычислительных ресурсов и некорректные настройки инструмента. Мы рассмотрим варианты того, как можно выяснить в чем проблема 18 2) Управление Задачами. Взаимодействие с задачами осуществляется либо с помощью главному меню, либо раскрывающегося меню

    19

    3) Просмотр задач. Любую задачу можно просмотреть и проверить лог-ошибок. Это позволит например выяснить какие настройки были использованы, а также просмотреть ошибки.

    20

    После того, как мы отправили задачку , она сначала появится в группе Running и будет там находится несколько минут, пока проводятся вычисления. 21 После завершения – задача будет перемещена в группу Complеted

    1. Получение и просмотр выходных файлов

    Теперь перейдем во вкладку файлового менеджера (FM). Нажмите на Files и выберите Balanophoraceae_project. В этой папке вы увидете созданные выходные файлы. По нашей задаче будет создано два файлы – оригинальный файл GenBank и csv файл. 23 Нам наиболее важна табличка csv . Просмотрим содержимое этого файла нажав на Preview кнопку в правом меню

    24

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

    25

    1. Интерпретация результатов и их анализ в R

    Теперь мы знаем какие гены и для каких видов были отскевенированы для Balanophoraceae и хотим выяснить, для каких генов отсеквенировано больше всего уникальных видов – ведь именно это позволят нам реконструировать дерево с наибольшим покрытием и запустить поиск BLAST. Для этого скопируйте функцию для R и проанализируйте файл Balanophoraceae.csv:

    gbSum <- function(in_csv_file){
        dd <- read.csv(in_csv_file, head = FALSE, stringsAsFactors = FALSE)
        rr <- data.frame(table(dd[, 1]))
        for (i in 1:dim(rr)[1]){
                rr[i, 3] <- dd[which(dd[, 1] == rr[i, 1]),][which(dd[which(dd[, 1] == rr[i, 1]),][, 3] == max(dd[which(dd[, 1] == rr[i, 1]),][, 3]))[1], 2]
                rr[i, 4] <- length(unique(dd[, 4][which(dd[, 1] == rr[i, 1])]))
            }
        colnames(rr)<-c("gene_name","frequency","reference_seq","unique_species")
        return(rr)
    }
    gbSum("Balanophoraceae.csv")

    Не забудьте указать директорию (setwd команда R)  в которую вы сохранили файлик Balanophoraceae.csv file 26   Из полученнного саммари мы видим следующее:
    1) три гена с максимальным видовым покрытием

    • 18S ribosomal RNA (15 уникальных видов)
    • 28S ribosomal RNA (7 уникальных видов)
    • 5.8S ribosomal RNA (7 уникальных видов)

    2) Референсные последовательности для данных генов: L24044.1 (18S) JN392881.1 (28S и 5.8S rRNA общая рефернсная последовательность)

    3) теперь мы должны приготовить текстовый файлик с разделителем табуляция где первая колонка – это референсная последовательность (относительно нее будет запущен BLAST) и название гена (для создания выходных файлов): 27

    Объявление и приглашение присоединиться к проекту InsideDNA

    Дорогие читатели блога!

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

    Мы разрабатывали новый сайт целых полтора года и теперь рады пригласить вас на борт нашего биоинформатического самолета! Встречайте – InsideDNA.me.

    logos

    Вот только небольшая часть функциональности сайта InsideDNA.me, к которой каждый пользователь получает доступ:

    • Доступ к более 600 наиболее популярным биоинформатическим инструментам (включая, например, TopHat, Bowtie2, OrthoMCL, samtools, bamtools, BEAST, phyml, abyss, SOAPdeNovo)
    • Все инструменты полностью настроены и готовы к работе через простой веб-интерфейс
    • 10 Гигабайт места на диске для хранения данных
    • Досуп к мощным серверам, включая сервера с более 200 Гб RAM и 32 ядрами (кроме того, число таких машин неограничено)
    • Возможность запросить подключение любой биоинформатической программы

    И многое другое, что можно будет увидеть, подключившись к нашему проекту. На данном этапе, мы можем подключить только 200 первых зарегистрировашихся пользователей.

    Теперь большую часть туториалов мы будем писать с использованием инструментов напрямую интегрированных в InsideDNA.me и многие вопросы и сложности, с которыми вы могли встечаться ранее при запуске скриптов и анализе своих данных, надеемся, будут решены быстро, удобно и помогут вам в написании научных статей и диссертаций.

    При регистрации мы просим новых пользователей заполнить анкету, чтобы понимать, кому оказывается полезным наш сайт. Итак – регистрируйтесь и пробуйте InsideDNA.me