BWA (англ.Burrows-Wheeler Aligner) — программный пакет для картирования коротких прочтений на большие референсные геномы (такие как, например, геном человека[2]), написанный китайским биоинформатиком Хенг Ли и англичанином Ричардом Дурбиным. Является одним из широкоиспользуемых алгоритмов выравнивания[3][4][5], а также рекомендуется для анализа данных производителями Illumina[6]. BWA состоит из трёх основных алгоритмов: BWA-BackTrack, BWA-SW и BWA-MEM. В основе алгоритмов BWA лежит преобразование Барроуза—Уилера, суффиксные массивы и алгоритм выравнивания Смита—Ватермана. Программный пакет умеет работать с длинными последовательностями на порядок быстрее, чем MAQ (англ.Mapping and Assembly with Qualities) при достижении аналогичной точности выравнивания[7].
Для образцов c платформ Illumina-SOLEXA, 454 Life Sciences, ионного полупроводникового секвенирования или секвенирования по Сэнгеру с длиной чтений («ридов» от англ.read — «чтение, прочтение»; отдельные фрагменты ДНК, полученные в результате работы секвенатора) больше или равных 70 пар оснований (здесь и далее — п.о.), BWA-MEM обычно предпочтительней. BWA-MEM и BWA-SW имеют схожие функции, такие как, например, поддержка длинных прочтений. BWA-MEM имеет лучшую производительность, чем BWA-BackTrack для чтений длиной 70—100 п.о. с платформы Illumina-SOLEXA.
Независимо от выбранного алгоритма программный пакет BWA сначала должен построить FM-индекс для референсного генома. Для создания индекса используется следующая команда: bwa index [-p prefix] [-a algoType] <in.db.fasta>. Далее к полученным данным применяются алгоритмы выравнивания, которые вызываются различными подкомандами: aln/samse/sampe для BWA-BackTrack, bwasw для BWA-SW и mem для BWA-MEM[8].
Хенг Ли, научный сотрудник Кембриджского института Броуда[9], начал писать первую часть кода 24 мая 2008 года и получил первоначальную стабильную версию уже 2 июня того же года. В этот же период профессор Tak-Wah Lam, первый автор статьи BWT-SW[10], сотрудничал с Пекинским геномным институтом в работе над SOAP2, преемником SOAP (пакет анализа коротких олигонуклеотидных последовательностей). SOAP2 вышел в ноябре 2008 года.
Алгоритм BWA-SW является новым компонентом BWA. Он был задуман в ноябре 2008 года и реализован спустя десять месяцев. По состоянию на 2015 год BWA-MEM (0.7.10) считался более предпочтительным в использовании для поиска снипов, чем bowtie (2-2.2.25) или Novoalign (3.02.12)[11].
Программный пакет BWA может быть использован в таких областях, как генотипирование для поиска однонуклеотидных полиморфизмов (снипов)[12]. Генотипирование позволяет проводить индивидуальное сравнение генотипов разных людей, с помощью него можно предсказывать предрасположенность людей к всевозможным заболеваниям, таким как, например, рак, диабет и многие другие. Так же пакет будет полезен в транскриптомике[13].
Начиная с версии 0.6.x, все алгоритмы могут работать с геномами, общая длина которых превышает 4 Гб[8]. Тем не менее, размер индивидуальной хромосомы не должен превышать 2 Гб. Каждый алгоритм вызывается отдельной командой, принимает входной файл через стандартный поток ввода («stdin») и возвращает результат через стандартный поток вывода («stdout»).
Сначала строится суффиксный массив, а затем генерируется BWT. Неопределённые нуклеотиды (обозначаются литерой «N») в референсном геноме преобразуются в случайные нуклеотиды, а в чтениях рассматриваются как несовпадения. Это может привести к регионам ложных совпадений («хитов») с неоднозначными основаниями. Но вероятность того, что это может произойти очень небольшая, учитывая относительно длинные чтения[14].
BWA-SW строит FM-индексы и для референсной последовательности и для чтений. Он неявно представляет собой референсную последовательность в виде префиксного дерева и чтения в виде детерминированного ациклического конечного автомата, который составляется из префиксного дерева для чтения с использованием динамического программирования. В отличие от BLAT и SSAHA2, BWA-SW находит сиды (от англ.seed — «семя»; точное совпадение нескольких пар оснований чтения[4]) в референсном геноме с помощью построения FM-индексов. Выравнивание расширяется тогда, когда сид имеет несколько вхождений в эталонной последовательности. Скорость достигается за счёт уменьшения ненужного расширения для высоко повторяющихся последовательностей. По сравнению со стандартным выравниванием Смита — Ватермана, временная сложность которого равна , BWA-SW имеет лучшую временную сложность так как он, как минимум, не медленнее, чем BWT-SW, сложность которого — [10]. В BWA-SW два выравнивания считаются различными в том случае, если длина зоны перекрытия по запросу составляет менее половины длины более короткого сегмента запроса.
Реализация BWA-SW принимает индекс BWA и FASTA запроса или FASTQ файла в качестве входных данных и выводит выравнивание в формате SAM. Входной файл, как правило, содержит много последовательностей («ридов»). Каждая последовательность обрабатывается, в свою очередь, с использованием нескольких потоков (если это возможно). Использование памяти при работе с данным алгоритмом составляет около 5 Гб для индексирования полных последовательностей генома человека[8]. Для коротких прочтений требуется приблизительно 3,2 Гб при использовании команды aln и 5,4 Гб при использовании sampe. Необходимый объём памяти для каждого запроса примерно пропорционален длине последовательности. Кроме того, BWA-SW способен обнаружить химеру, потенциально вызванную структурными изменениями или эталонными ошибками сборки ридов, которые могут представлять угрозу для BLAT и SSAHA2. BWA-SW отличается от BWT-SW в нескольких аспектах[15]: BWT-SW гарантирует, что найдёт все локальные совпадения, в то время как BWA-SW представляет собой эвристический алгоритм, который может пропустить истинные совпадения, но работает гораздо быстрее.
BWA-MEM является последним реализованным алгоритмом пакета BWA. Он использует стратегию, в которой алгоритм выравнивания прочтений выбирается автоматически в каждом конкретном случае: локальное или глобальное выравнивание. Поддерживаются парноконцевые прочтения и выравнивания химерных прочтений. Алгоритм устойчив к ошибкам секвенирования и применим к широкому диапазону последовательности длиной от 70 п.о. до нескольких мегабаз. А для картирования последовательностей длинной 100 п.о. BWA-MEM показывает лучшую производительность, чем некоторые передовые алгоритмы выравнивания прочтений[16].
Для каждого выравнивания программный пакет BWA вычисляет оценку качества картирования — Phred масштабированную вероятность того, что выравнивание является некорректным[8] (англ.Phred quality score — оценка качества прочитанного нуклеотида ДНК). Необходимость статистической оценки картирования возникает из-за того, что риды получаются маленькие (40—100 п.о), а геномы достаточно большие, поэтому выравниваний может быть несколько. Тем не менее, необходимо понять, какое из них наиболее информативное. И по сути Phred-оценка представляет собой вероятность того, что рид выравнен неправильно. Вычисляется она по формуле , где q — качество прочтения рида[17]. Алгоритм аналогичен оценке для MAQ, за исключением того, что в BWA предполагается, что истинное совпадение всегда может быть найдено. Такое изменение было сделано в связи с тем, что формула MAQ переоценивает вероятность пропуска истинного совпадения, что приводит к заниженному качеству картирования. Моделирование показывает, что BWA может переоценивать качество картирования из-за этой модификации, но отклонение относительно невелико.
Пример сценария использования программного пакета BWA от Голландского вычислительного центра выглядит следующим образом[18]:
Шаг 1 — индексирование генома (∼3 CPU часа для генома человека)
bwa index -a bwtsw ref.fa
,где ref.fa — название файла с референсным геномом.
Шаг 2а — генерирование выравнивания в координатах суффиксного массива
bwa aln ref.fa read1.fq.gz > read1.sai
bwa aln ref.fa read2.fq.gz > read2.sai
,где read1.fq.gz — файл с прямыми чтениями, read2.fq.gz — обратными чтениями.
Применяем параметр -q15, если качество меньше на 3’-концах чтений.
Шаг 3а — генерирование выравнивания в SAM формате для парноконцевых чтений
bwa sampe ref.fa read1.sai read2.fq.gz > aln.sam
Шаг 4a — генерирование выравнивания в SAM формате для непарноконцевых чтений
bwa samse ref.fa read1.sai read1.fq.gz
Шаг 2б — BWA-SW для длинных чтений
bwa bwasw ref.fa long-read.fq.gz > aln-long.sam
,где long-read.fq.gz — файл с длинными чтениями.
Для дальнейшей совместимости полученного выравнивания со стандартными программами для анализа NGS данных, таких как, например, Pikard Tools или GATK, при вызове программы следует указывать так называемую группу прочтения (англ.read group — набор данных, полученных с помощью одного режима работы секвенсора) добавлением флага -r[19]. Сама группа прочтения для образца с платформы Illumina-SOLEXA выглядит примерно следующим образом: -r ‘@RG\tID:1\tSM:S000336\tPL:ILLUMINA\tLB:L00000336\t PU:C3A7NACXX:1’[19], где ID — это уникальный идентификатор данного образца, SM — это названия образца, Pl — это платформа для секвенирования, LB — название библиотеки, а PU — системный идентификатор[20].
На практике BWA является основным инструментом, предназначенным для выравнивания коротких ридов[21][22][23][24]. В настоящее время BWA применяется в следующих исследованиях:
Секвенирование и сборка генома.
Поиск однонуклеотидных полиморфизмов.
Поиск делеций и вставок.
Основными преимуществами использования программного комплекса BWA являются следующие его достоинства[7][25]:
↑BWA Usage(англ.). Голландский вычислительный центр. Архивировано из оригинала 5 апреля 2017 года.
↑ 12Bukowski Robert, Sun Qi, Wang Minghui.Variant calling. — Bioinformatics Facility, Institute of Biotechnology. — С. 22. Архивировано 6 апреля 2017 года.
, BWA-SW имеет лучшую временную сложность так как он, как минимум, не медленнее, чем BWT-SW, сложность которого — 
[10]. В BWA-SW два выравнивания считаются различными в том случае, если длина зоны перекрытия по запросу составляет менее половины длины более короткого сегмента запроса.
, где q — качество прочтения рида[17]. Алгоритм аналогичен оценке для MAQ, за исключением того, что в BWA предполагается, что истинное совпадение всегда может быть найдено. Такое изменение было сделано в связи с тем, что формула MAQ переоценивает вероятность пропуска истинного совпадения, что приводит к заниженному качеству картирования. Моделирование показывает, что BWA может переоценивать качество картирования из-за этой модификации, но отклонение относительно невелико.