суббота, 11 июня 2016 г.

Парадигма картирования seed-and-vote: итоговые файлы в формате SAM.

Успешная работа такого элайнера как Subjunc, равно как и любого другого элайнера, завершается генерацией итогового файла, содержащего информацию о результатах картирования ридов. Исходным форматом такого файла является формат SAM (акроним от англ. Sequence Alignment/Map), но для удобства записи, хранения и дальнейшего применения конечному пользователю обычно предлагается его бинарный (двоичный) BGZF-сжатый вариант BAM (акроним от англ. Binary Alignment/Map). В данном сообщении мы рассмотрим ключевые спецификации формата SAM.
Формат SAM был разработан международной группой ученых из Wellcome Trust Sanger Institute (Cambridge, UK), Broad Institute of MIT and Harvard (Cambridge, USA), Beijing Institute of Genomics (Beijing, China), Department of Computer Science at University of California (Los Angeles, USA), Department of Biology at Boston College (Chestnut Hill, USA) и Center for Statistical Genetics at University of Michigan (Ann Arbor, USA), а также членами Data Processing Subgroup проекта 1000 Genomes. Публично формат был анонсирован и описан в 2009 году в специализированной статье “The sequence alignment/map format and SAMtools, опубликованной в 16-ом номере 25-го тома журнала Bioinformatics. Там же был описан пакет SAMtools, ныне являющийся одним из наиболее популярных свободно распространяемых пакетов программ для работы с SAM/BAM файлами.
Формат SAM является разновидностью текстового формата хранения данных, в котором отдельные поля данных разделены с помощью табуляции. Все данные в файлах такого формата собраны в двух разделах: “Заголовок” (header section) и раздел “Картирование” (alignment section). Первый из этих разделов не является обязательным и может отсутствовать, второй же является обязательным.

Описание формата SAM.
Раздел “Заголовок”.
Каждая строка в этом разделе начинается с символа @, за которым идет двухбуквенное кодовое обозначение группы данных, записанных в строке. Данные в каждой строке разделены табуляцией и записаны в формате “ярлык:значение” (“TAG:VALUE”). Для ярлыков применяется двухбуквенная кодировка, причем для общепринятых и устоявшихся ярлыков используются прописные буквы, а для пользовательских (предложенных самим пользователем) – строчные. Значения ярлыков могут содержать различную информацию и по разному записаны, что во многом определяется групповой принадлежностью данных. Ниже дана расшифровка основных групп данных и кодировок, принятых в заголовке файла формата SAM.

Ярлык
группы данных
Ярлык
данных
Значение ярлыка
и его описание
@HD

Строка, посвященная заглавным данным.
VN*
Версия формата. Допускаемые значения: /^[0-9]+\.[0-9]+$/
SO
Способ сортировки картированных ридов. Допускаемые значения: unknown (по умолчанию), unsortedqueryname и coordinate. При сортировке по координатам в качестве основного ключа используется поле RNAME (описание этого поля, а также ряда других полей, упомянутых далее, дано в разделе “Картирование”), а порядок сортировки прописывается в заглавной строке @SQ. Второстепенным ключом для сортировки является поле POS. При равности полей RNAME и POS риды сортируются произвольно. Те риды, которые в поле RNAME помечены символом * ставятся после ридов, помеченных любым другим символом.
GO
Способ группировки картированных ридов. Допускаемые значения: none (по умолчанию), query (похожие риды сгруппированы по значению поля QNAME) и reference (похожие риды сгруппированы по значению полей RNAME и POS).
@SQ

Строка, посвященная референс-последовательности, относительно которой проводилось картирование ридов. Поскольку референс-последовательностей может быть много (например, все хромосомы человека), то и строк @SQ также может быть много, причем очередность этих строк определяет порядок сортировки ридов в файле формата SAM.
SN*
Название референс-последовательности. Если строк @SQ много, то каждая из них должна иметь уникальное значение ярлыка SN.
LN
Длина референс-последовательности, которая должна лежать в диапазоне от 1 до 231 – 1.
AS
Идентификатор сборки генома.
M5
Контрольная сумма референс-последовательности, записанная прописными буквами без пробелов.
SP
Видовая принадлежность референс-последовательности.
UR
Унифицированный идентификатор ресурса, где содержится референс-последовательность. Обычно начинается с метки одного из стандартных протоколов обмена данными http или ftp, а в отсутствии таковой принимается, что идентификатор указывает путь к референс-последовательности в файловой системе локального компьютера.
@RG

Строка, посвященная группе ридов. Поскольку групп ридов может быть много, то и строк @RG также может быть много, причем очередность этих строк может быть произвольной.
ID*
Уникальный идентификатор группы ридов. Если строк @RG много, то каждая из них должна иметь уникальное значение ярлыка ID. При объединении нескольких SAM файлов, во избежание совпадений, идентификаторам групп присваиваются новые “сквозные” уникальные значения.
CN
Название учреждения, где было осуществлено секвенирование.
DS
Описание группы ридов.
DT
Дата, когда было проведено секвенирование (указанная в соответствии со стандартом ISO8601 или в произвольной форме типа дата/время).
FO
Последовательность потоков нуклеотидов, запускаемых при чтении каждого из ридов. Потоки имеют символьную запись в формате /\*|[ACMGRSVTWYHKDBN]+/.
KS
Массив нуклеотидов, использованных для прочтения ридов.
LB
Библиотека.
PG
Программы, использованные для обработки группы ридов.
PI
Прогнозируемая медиана размера вставки.
PL
Платформа/технология, использованная для секвенирования. Допускаемые значения: CAPILLARY, LS454, ILLUMINA, SOLID, HELICOS, IONTORRENT, ONT, и PACBIO.
PM
Произвольный текст, содержащий более детальное описание платформы/технологии, использованной для секвенирования.
PU
Уникальный идентификатор загрузочного модуля.
SM
Название образца, предназначенного для секвенирования.
@PG

Строка, посвященная программе, формирующей запись в SAM файле. Если запись в SAM файле осуществлялась несколькими программами, то строк @PG может быть много.
ID*
Уникальный идентификатор программы.
PN
Название программы.
CL
Командная строка, использованная для запуска и контроля за работой программы.
PP
Уникальный идентификатор предыдущей программы, работавшей над записью в SAM файле. При последовательной работе множества программ данный ярлык создает их упорядоченный список, в котором первой стоит программа, осуществившая последнюю запись в SAM файл.
DS
Описание программы.
VN
Версия программы.
@CO

Строка, содержащая произвольный текстовый комментарий. Допускается неупорядоченная запись множества таких строк. Причем это единственная группа данных в разделе “Заголовок”, которая не разбивается табуляцией на отдельные поля.

Если есть соответствующая заглавная строка данных, то те ярлыки, которые в выше представленном описании помечены символом *, являются обязательными. Остальные ярлыки не являются обязательными и могут быть пропущены. Кроме того, пользователь может добавлять в раздел “Заголовок” любые другие ярлыки, поля и значения. В этом случае рекомендуется лишь использовать не прописные, а строчные буквы, что бы различать пользовательские обозначения от обозначений, закрепленных в официальной спецификации формата SAM.
Ниже в качестве примера приведен небольшой фрагмент раздела “Заголовок” SAM файла, сгенерированного элайнером Subjunc:

@HD
VN:1.0
SO:coordinate

@SQ
SN:chr1
LN:248956422

@SQ
SN:chr10
LN:133797422

@SQ
SN:chr11
LN:135086622

@PG
ID:subread
PN:subread
VN:Rsubread 1.20.2

Раздел “Картирование”.
Этот раздел состоит из множества строк, по одной на одно картирование. При этом одному риду может быть отведено несколько строк (например, если он картирован по множеству мест в референс-последовательности или же имеет химерную природу). Каждая строка разбита на 11 обязательных полей, идущих в том порядке, в котором они описаны ниже.

Поле
Формат значения
поля (регулярное выражение и/или диапазон)
Описание поля
QNAME
[!-?A-~]{1, 254}
Название матрицы (query template name), которая была полностью или частично прочитана во время секвенирования. Понятие матрица в данном случае тождественно понятию фрагмент ДНК, но не всегда тождественно понятию рид, так как во время секвенирования с одного и того же фрагмента ДНК может быть считано несколько ридов (например, в случае с парными ридами). Если так, то строки ридов, полученных с одного и того же фрагмента ДНК, будут иметь одно и то же значение QNAME.
FLAG
[0, 216-1]
Битовая кодировка свойств проведенного картирования. Каждый SAM FLAG соответствует одному или нескольким свойствам. Так, SAM FLAG, равный 4, обозначает, что рид не картирован, а SAM FLAG, равный 65, свидетельствует о том, что рид принадлежит паре и что в этой паре он является первым. Интерпретация основных SAM FLAGs дана в официальной спецификации формата SAM, доступной через веб-сервис GitHub на странице проекта SAMtools. Интерпретацию же любых побитовых комбинаций SAM FLAGs можно получить через ряд специализированных веб-сервисов, например, Picard.
RNAME
\*|[!-()+-<>-~][!-~]*
Название референс-последовательности, по которой осуществлено картирование. В том случае, когда в разделе “Заголовок” есть одна или несколько строк @SQ, то хотя бы один ярлык SN должен иметь такое же название. Если поле RNAME имеет значение *, то поля POS и CIGAR (см. далее) не могут быть заполнены актуальными значениями.
POS
[0, 231-1]
Координата (в референс-последовательности) крайнего левого успешно картированного нуклеотида рида. В самой же референс-последовательности первый нуклеотид имеет координату 1. Если рид не картирован, то полю POS присваивается значение 0, а поля RNAME и CIGAR принимают значение *.
MAPQ
[0, 28-1]
Качество картирования рида, выраженное в баллах (MAPQ scores). Балы рассчитываются по формуле –10 * log10(вероятность ошибочного картирования). Так, если вероятность корректного картирования равна 0.99, то MAPQ score будет равен –10 * log10(1 – 0.99) = 20. Поскольку при таком подходе могут получаться дробные значения, то принято их округлять до ближайшего целочисленного значения. Очевидно, что предельное значение поля MAPQ задается желаемой максимальной точностью оценки качества картирования. При этом у разных элайнеров используются разные предельные значения поля MAPQ. Так, например, предел, используемый BWA, равен 37, Bowtie 2 – 42, а Subjunc – 59. Максимум, который допустим по спецификации формата SAM, равен 255, но в действительности это значение присваивается полю MAPQ только тогда, когда по тем или иным причинам реальное значение качества картирования не может быть получено.
CIGAR
\*|([0-9]+[MIDNSHPX=])+
Краткая запись результатов выравнивания рида относительно референс-последовательности. Приняты следующие обозначения: M – нуклеотид совпадает, I – нуклеотид вставлен, D – нуклеотид удален, N – нуклеотид (в референс-последовательности) “перепрыгнут” (в тех случаях, когда относительно референс-последовательности выравнивается рид, перекрывающий в молекуле РНК экзон-экзонный переход), S – рид содержит “мягкий клипс” (последовательность с 5’- или 3’-конца рида, которая “свисает”, так как нигде не картируется, и может быть обрезана), H – рид содержит “жесткий клипс” (последовательность с 5’- или 3’-конца рида, которая “свисает”, но не может быть обрезана, так как картируется по другому участку референс-последовательности в другой ориентации), P – в рид добавлен ”пустой” нуклеотид *, благодаря чему рид был успешно картирован. С учетом описанных обозначений CIGAR запись вида, например, 20M2I79M, обозначает, что у рида полностью совпадают с референс-последовательностью первые 20 нуклеотидов, затем идет небольшая 2-нуклеотидная вставка, после которой находится еще 79 полностью совпадающих с референс-последовательностью нуклеотидов.
RNEXT
\*|=|[!-()+-<>-~][!-~]*
Название референс-последовательности, по которой был картирован следующий рид. Под следующим ридом подразумевается рид, который идет сразу же за текущим и который был считан с того же фрагмента ДНК (или матрицы, см. поле QNAME), что и текущий рид. Если текущий рид является последним из всех ридов, считанных с данного фрагмента ДНК, то следующим ридом для него будет первый рид из всех ридов, считанных с данного фрагмента ДНК (например, в случае парных ридов для первого рида следующим ридом является второй рид пары, а для второго рида следующим будет первый рид пары). В том случае, когда в разделе “Заголовок” есть одна или несколько строк @SQ, то хотя бы один ярлык SN должен иметь такое же название. Это поле может принимать также значение * если необходимая информация недоступна, или значение = если RNEXT идентично RNAME.
PNEXT
[0, 231-1]
Первая координата (в референс-последовательности) следующего рида. Эта координата равна значению поля POS самого следующего рида. В том случае, когда необходимая информация недоступна, полю PNEXT присваивается значение 0.
TLEN
[-231+1, 231-1]
Длина фрагмента ДНК (или матрицы, см. поле QNAME), с которого был считан данный рид. Этому полю присваивается значение 0 при условии что: 1) поле RNEXT имеет значение *; 2) с фрагмента ДНК было считано два (или более) рида, но они картировались по разным референс-последовательностям; 3) с фрагмента ДНК был считан один рид, но во время выравнивания он был разбит на несколько сегментов и эти сегменты были картированы по разным референс-последовательностям. Если же с данного фрагмента ДНК было считано два (или более) рида (или же один рид, но он во время выравнивания был разбит на несколько сегментов) и все они картировались по одной и той же референс-последовательности, то за длину фрагмента ДНК принимается количество нуклеотидов от крайнего левого до крайнего правого из картированных нуклеотидов. При этом значение поля TLEN будет положительным, если текущий рид является крайним левым из всех ридов, считанных с фрагмента ДНК, и отрицательным (но равным по модулю), если речь идет о крайнем правом риде из всех ридов, считанных с фрагмента ДНК.
SEQ
\*|[A-Za-z=.]+
Нуклеотидная последовательность рида. Если нуклеотидная последовательность рида не сохранена, то этому полю присваивается значение *.
QUAL
[!-~]+
Качество прочтения нуклеотидов в риде, записанное символами ASCII. Длина символьной строки равна длине рида из поля SEQ. Значение для поля QUAL берется из FASTQ файла. Если эта информация не сохранена, то полю присваивается значение *.

Помимо обязательных полей каждая строка может содержать ряд дополнительных (но не обязательных) полей. Например, поле OQ используется для хранения символьной строки с исходным качеством прочтения нуклеотидов в риде, если проводилась рекалибровка качества. Подробное описание всех дополнительных полей раздела “Картирование” можно найти в официальной спецификации формата SAM, доступной через веб-сервис GitHub на странице проекта SAMtools.
Ниже в качестве примера приведен небольшой фрагмент раздела “Картирование” SAM файла, сгенерированного элайнером Subjunc:

SRR1145838.42031698
99
chr1
14656
59
100M1S
=
14693
138
CCCAAGTCTG
CCCFFFFFHH
HI:i:1
NH:i:1
NM:i:0
SRR1145838.27603184
163
chr1
14669
57
8M1D93M
=
14718
149
CGGGGGGAAA
CCCFFFDDDD
HI:i:1
NH:i:1
NM:i:3
SRR1145838.35898901
163
chr1
14669
58
101M
=
14755
187
CTGGGGGGGA
@C@DDFFCC3
HI:i:1
NH:i:1
NM:i:1
SRR1145838.15702169
89
chr1
14669
58
101M
*
0
0
CTGGGGGGAA
?DDDDDDCDD
HI:i:1
NH:i:1
NM:i:1
SRR1145838.2006935
145
chr1
16146
59
101M
chr16
15765
0
ATCCTGTGTGC
:@DDBACDDD
HI:i:1
NH:i:1
NM:i:0
SRR1145838.8452425
147
chr1
16239
59
1S72M296N28M
=
16182
-157
CCACACGAGCC
(88A0@DDDD
HI:i:1
NH:i:1
XS:A:-

Ради компактности в этом примере сиквенсы ридов и соответствующие поля качества прочтения нуклеотидов сокращены до 10 символов.

Комментариев нет:

Отправить комментарий