Практическая биоинформатика

    Обнаружил жесткую нехватку информации по биоинформатике в русском сегменте. Не знаю, уж, востребована она или нет, но хочу предоставить на суд читателя вводную часть, которую можно назвать практическая биоинформатика, которой мне очень не хватало для ознакомления с предметом. В этой главе я хочу описать путь, который пришлось пройти мне до настоящего момента, когда я уже не шарахаюсь от фраз: вот вам FASTQ файл и постройте мне bed graph для genome browser. Чтобы в дальнейшем вести разговор об интересном, хочу по диагонали пройтись по определениям и программам первичной обработки данных, без которых трудно говорить на одном языке.

    Для начала немного определений. Будем считать, что хромосомы это одномерные координатные оси, от единицы до, примерно, 10е8. Длина оси зависит от длины хромосомы. Каждая точка оси, это целое число. Биологи и химики провели большое количество экспериментов и благодаря им смогли уже с большой точностью (около 90%) описать части хромосом. Эти описания называются аннотацией. Аннотация содержит информацию о длине хромосом, о координатах отдельных участков хромосом, наиболее известны из них ген, intron, exon. Этих участков огромное количество, но основное их свойство заключается в том, что это отрезки или набор отрезков расположены на координатной оси. Одни отрезки могут включать в себя другие или как-либо пересекаться. Вот набор сайтов на котором можно посмотреть аннотации геномов человека, мышки и тому подобных.

    Биологи с химиками проводят свои эксперименты и в результате их операций над клетками получают раствор, в котором содержатся относительно небольшие куски ДНК или РНК (не очень хочется вдаваться в подробности разницы или одинаковости, просто последовательность нуклеотидов). Этот раствор пропускают через аппаратуру для секвенирования, выход которой — небольшие строки. Эти строки являются концами кусков ДНК или РНК, находившимися в растворе. Длина полученных из аппаратуры строк составляет всего 36-50 bases (длина строки в нуклеотидах) иногда больше, но на текущий момент вроде не более 200. Эти отрезки, полученные из секвенирующей аппаратуры и определенные, последовательностью нуклеотидов, называются ридами (от английского reads — “считывание”). Стоит отметить, что риды характеризуются только последовательностью нуклеотидов, а не расположением на геноме. Иногда эти последовательности дополнены строкой вероятностей, ставящей в соответствие позиции нуклеотида его вероятность нахождения на этой позиции. FASTA файл — это файл без вероятностей, FASTQ — это файл с вероятностями.

    Далее, в зависимости от того, что стало результатом эксперимента — куски ДНК или РНК, проводят одну из двух методик секвенирования ChIP-seq или RNA-seq, соответственно. Более подробно о них рассказано тут http://en.wikipedia.org/wiki/DNA_sequencing.

     После того, как дорогие машины по секвенированию отработают и выдадут результат в FAST A/Q файл, необходимо найти получившиеся последовательности в геноме. Для ChIP-seq очень быстрый поиск, возможный в домашних условиях, делает программа bowtie, находила миллионы ридов всего за 5 минут. Т.е. она ищет вхождение строки длиною 36-50 символов, состоящей из, как минимум, четырех буквенного алфавита в строках общей длиною до 10е9. Почему использован оборот как минимум: помимо стандартного использования A/T G/C алфавита часто добавляют символ N, заменяющий любую возможную букву (более детально http://genome.ucsc.edu/FAQ/FAQdownloads.html#download5). У программы множество параметров. Так, например, можно разрешить в строке(риде) ошибку или две (максимум до 3х ошибок). Она может искать не одно вхождение рида в геном, а много. Может хитро сортировать данные, например, если рид с одной ошибкой нашелся однозначно в геноме, а с двумя или тремя ошибками нашелся во многих частях генома, то вывести только первый результат. Этот процесс нахождения ридов в геноме называется мэппингом от английского mapping (еще его называют alignment). Алгоритм такого быстрого поиска очень интересен, но этому можно посвятить отдельную статью или найти ссылку на английскую статью, с сайта разработчика, где рассказывается как алгоритм bzip2 натолкнул их на такое решение. Программ по мэппингу и сайтов, где его делают онлайн, много; по ключевым словам blast, eland +genome в гугле можно найти дополнительную информацию.

    Для RNA-seq процедура немного сложнее, для него выполняют ChIP-seq мэппинг, а далее обрабатывают те риды, которые не нашли себя при ChIP-seq. Хорошая программа, которая активно использует bowtie для своей работы, называется tophat. В результате сплайсинга и порожденного им изоформизма части рида могут находиться в разных местах генома. Например, первые 15 символов могут приходиться на один участок генома, а другие 11 — на другой. Это разделение рида на части, приходящееся на конец одного exon и на начало другого, называют splice-junctions. Tophat позволяет их находить, а также определяет новые возможные изоформы существующих генов.

    Результатом работы этих программ является sam/bam файл, в котором содержится информация из FAST A/Q файла плюс информация о координатах на соответствующих спиралях хромосомы. Процесс от библиотеки до sam/bam файла, часто называют pipeline процедурой и во многих лабораториях стоит на потоке, поэтому стоит поинтересоваться, что за параметры и версии программ установлены по умолчанию. В общем, это то, на чем заканчивается вступительная часть, далее наступает очередь исследований. Можно отметить, что с этого момента у нас есть данные, находящиеся с некоторой достоверностью в одних и тех же координатах: координатная ось, аннотация, риды с соотвествующими координатами.

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

    Если вас заинтересовала вводная часть, то на каждом этапе я могу остановиться подробнее. К сожалению, простая вводная часть уже занимает несколько страниц сухого текста без картинок, поэтому описывать написанные программы и механизмы в этой главе не решился. Меня самого больше всего привлекает последний абзац, где я мельком упомянул про статистику. Хочется освятить существующие библиотеки и механизмы для работы с такими данными. Сюда можно приобщить методы, входящие в datamining (разные типы кластеринга), которые описаны на хабре. Как применить пуассоновское распределение для анализа данных без контроля, как применить сложную цепочку пуассона, ф-теста, для нахождения участков обогащенных ридами на геноме (Diarac delta function)? Полезны ли готовые библиотеки для работы с интервалами boost.intervals, boost.icl?

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

    На текущем этапе своей работы мы пытались найти параметры, по которым можно фильтровать риды, как на интересные для исследования, так и на шумовые. Задача довольно нетривиальная при отсутствии контроля. В дальнейшем решено в библиотеку для секвенирования добавлять контроль, который позволит измерять уровень ошибки, но статистика все равно не отпустит.
Share post
AdBlock has stolen the banner, but banners are not teeth — they will be back

More
Ads

Comments 15

    +4
    Эта тема крайне интересна, и на русском материалов я не особо видел. Был бы крайне рад продолжению.
      +1
      Я только начал работать в сфере биоинформатики, пока пол года. Стараюсь собрать все мысли воедино и оформить в статьи. По молекулярной биологии брал классы, сейчас хожу на классы functional genomic. Часть вещей подчерпнул из статей. Есть огромное количество инетересных статей описывающие математические подходы для решения ряда задач связанных с анализом ридов после сиквенирования. Вот разрваюсь толи статьи осветить, толи дать примеры на С++ как с sam/bam файлами работать. Ибо в дальнейшем всеравно придется с данными работать.
      Надо собраться и написать и те статьи и те.
      Меня ещё гложет вопрос, хочется привлечь математиков и программистов к обсуждению. А тема подходящая биотехнологии. Не нашел на хабре Биоинформатика.
        0
        Пождите, с примерами С++. Не скажу за других, но мне бы для начала было интересно, да не только интересно, но и необходимо для понимания материала, узнать доступным языком такие вопросы:

        — что есть ДНК?
        — можно ли сказать, что ДНК несет бинарный код, подобно компьютерам? Пары нуклеотидов AT-CG подобны 0 и 1 у компьютера?
        — что есть хромосома? Скрученная ДНК?
        — Что есть ген? Фрагмент ДНК имеющий некую определенную смысловую нагрузку?

        Ну в общем если напрячься я таких базовых вопросов могу еще поднасыпать. :)
      0
      больше красивых картинок!!!
        0
        Такой вот был пост. Там ссылки на лекции (в т.ч. видео).
          0
            0
            Видео по теме
            Обзорная лекция по биоинформатике. Павел Певззнер.

            Директор Центра вычислительной масс-спектрометрии Национального института здоровья. Профессор факультета компьютерных наук и инженерии Университета Калифорнии (Сан-Диего). Победитель конкурса президентских мегагрантов в РФ (2010)

            Прочитана для старшеклассников.
            www.lektorium.tv/lecture/?id=13384
              +1
              Имею практический вопрос

              Есть ли какой-нибудь публичный API для доступа к геномным базам данных, того же пабмеда? Ваяю сейчас на .met одну геномную тулзу и прикидываю на будущее. Идеально конечно чтобы API было прям на готовые методы, типа BLAST.
                0
                .net :)
                  0
                  Я, к сожалению, пока к pubmed ничего не писал. Работаю в основном с UCSC genome browser. Про него знаю, но не много. Но про API ничего не слышал всеравно, просто SQL к BD.
                    0
                    Ну sql это тоже api в каком-то смысле :)
                      0
                      Согласен, но в моем случае, я локально установил базу. А под API наверное подразумевается, возможность обращаться к родным ресурсам.
                        0
                        >А под API наверное подразумевается, возможность обращаться к родным ресурсам.

                        Ага.

                        Хотя бы через SQL
                    +1
                    0
                    Спасибо, интересно очень. Давно хочу разобраться в этой теме, но как то тяжко построить в голове четкую картину. Нужно больше почитать и что-то попробовать. На topcoder сейчас много соревнований с призами на тему биоинформатики.

                    Only users with full accounts can post comments. Log in, please.