СБОРКА ГЕНОМА DE NOVO

ПОДГОТОВКА ЧТЕНИЙ

В данном практикуме работа ведется с проектом по секвенированию бактерии Buchnera aphidicola (код доступа SRR4240389). Это короткие (длины 36) чтения, полученные по технологии Illumina. Файл был скачен из базы данных ENA, рапспакован с помощью команды:

gunzip SRR4240389.fastq.gz

1. В первую очредь из полученного файла были удалены адапторы с помощью программы trimmomatic. Для этого использовалась следующая команда (где adapters.fasta - файл, содержащий последовательности адаптеров для Illumina):

java -jar /usr/share/java/trimmomatic.jar SE -phred33 SRR4240389.fastq SRR4240389_clean.fastq ILLUMINACLIP:adapters.fasta:2:7:7

После работы этой программы из 12950609 ридов, изначально содержавшихся в файле осталось 12947338 (99,97%), остальные 3271 (0,03%) были удалены, при этом по объему файл сократился на 1 Мбайт.

2. После этого была проведена операция отрезания от ридов концевых нуклеотидов с Quality score меньше 20 (TRAILING:20) и длиной менее 30 (MINLEN:30) с помощью той же программы trimmomatic, команды:

java -jar /usr/share/java/trimmomatic.jar SE -phred33 SRR4240389_clean.fastq SRR4240389_cleaner.fastq TRAILING:20 MINLEN:30

Как было показано ранее, такая операция существенно увеличивает качество прочтений. В результате работы программы из оставшихся 12947338 ридов было удалено более половины 6603310 (51,00%), осталось 6344028 (49,00%) прочтений. Размер файла сократился с 1386 Мбайт до 639 Мбайт.

ПОДГОТОВКА K-МЕРОВ

Для подготовки k-меров (или слов длины k) длины 29 и 25 использовалась программа velveth. На выходе ее работы появилась новая папка с тремя файлами: log, roadmaps, sequences.

velveth velveth 29 -fastq -short SRR4240389_cleaner.fastq
velveth velveth1 25 -fastq -short SRR4240389_cleaner.fastq

CБОРКА НА ОСНОВЕ K-МЕРОВ

С помощтю программы velvetg из полученных файлов с k-мерами были осущетвлены построение ориентированного графа де Брайна и сборка генома:

velvetg velveth
velvetg velveth1

Програма дает на выход файлы Graph, Pregraph, Lastgraph, contigs.fa (содержит последовательности всех контигов), stats.txt (содержит информацию о всех контигах, состоящих на вершинах графа). В Таблице 1 приведена основная информация об этих двух сборках:

Таблица 1. Характеристика сборок.

Длина К-мера	N50	Число контигов	Максимальная длина контига	Длина сборки
29	2960	1150	15093	654172
25	2563	4559	17211	704335

Как видно из Таблицы 1. N50 на 400 больше у сборки с К-мером длины 29, что идет вразрез с длиной контига, большей для сборки с К=25, число вершин графа также больше для 25, ровно как и число ридов. Далее был произведен более подробный анализ контигов: в Таблице 2 приведены средние значения для контигов, посчитанные с помощью Excel и информация о самых длинных контигах.

Таблица 2. Характеристика контигов (средние значения и самые длинные).

29			25
ID	Длина	Покрытие	ID	Длина	Покрытие
Среднее значение	568,85	598,21	Среднее значение	154,49	363,46
Медиана	29	20,038	Медиана	6	8,33
36	15093	37,135096	19	17211	43,595143
2	11042	47,141188	62	10842	59,173492
62	10741	38,155386	79	9315	56,203865

В Таблице 3 приведены данные об аномальном покрытии встречающимся в обеих сборках (превышение среднего значения в пять раз, или значение в пять раз меньше среднего), а также о контигах с аномально большим покрытием. Данные были взяты из файла stats.txt папки velveth, образованного после сборки k-меров. Как видно, существует множество контигов с аномально низким покрытием. Контиги с аномально высоким покрытием имеют длину 1, и даже не содержатся в файле contigs.fa, однако информация о них оказалась вместе с информацией о более крупных контигах в файле stats.txt.

Таблица 3. Аномальные покрытия, информация о двух континах с аномально высоким покрытием.

Сборка	Число контигов с аномально низким покрытием	Число контигов с аномально высоким поккрытием
29	1036	2
25	3408	49

Номер контига	Длина	Покрытие
1	1	597662
696	1	6436

АНАЛИЗ

Далее с помощью программы megablast каждый из трёх самых длинных контигов из Таблицы 2 был сравнен с хромосомой Buchnera aphidicola (GenBank/EMBL AC — CP009253). В таблице 4 приведены соответствующие результаты.

Таблица 4. Сравнение сборки k-29 c референсным геномом.

ID контига	Координаты на хромосоме	Число Гэпов	% Идентичности	Max score	Total score	Query cover	E-value
2	534945- 546014	285/11212(2%)	8834/11212(79%)	7267	7267	99%	0.0
36	621055- 620926	246/6173(3%)	4681/6173(76%)	2907	5147	63%	0.0
62	584329- 587055	110/2778(3%)	2101/2778(76%)	1279	2378	49%	0.0

В первом случае было построено лишь одно выравнивание, для 62, 36 контига было построено 2 выравнивания, выбрано для таблицы было с большим max score. Ниже приведены иллюстрации полученных выравниваний, на которых видно, что в целом контиги ровно ложаться на хромосому, присутствую также небольшие делеции.

Рис.1 Выравнивание контига 2	2 контиг на 99% покрывается хромосомой, совпадающие фрагменты практически идентичны, без каких-либо перестроек (без инфверсий, крупных делеций и т.п).
Рис.2 Выравнивание контига 36	Контиг 36 покрывается хромосомой лишь на 63%, из Рис. 2 видно, что аналогично первому случаю, хромосома ложится на контиг ровно, однако имеет место крупная делеция. Второе выравнивание покрывает хромосому полностью, без крупных делеций, однако Score для того выравнивания ниже, %Identities фактически идентичны (2239; 78% против 2907; 76%). Координаты выравнивания на хромосоме для двух выравниваний также отличаются своими границами (621055 to 627104 для представленного в таблице, 617419 to 620926 для не представленного).
Рис.3 Выравнивание контига 62	Для контига 62 видно, что он иныертирован относительно хромосомы, присутствует довольно крупная делеция, которой можно частично объяснить низкий процент покрытия контига хромосомой. Второе построенное выравнивание ниже по Score и %Identities приведенного в таблице (1099 и 75% против 1279 и 76%). Отличаются и кординаты выравнивания на хромосоме (584329 to 587055 и 580083 to 582686).

СБОРКА ПОЛОВИНЫ РИДОВ

С помощью программы tail была удалена первая половина очищенного от остатков адапторов, от концевых низкокачественных нуклеотидов и от коротких ридов файла SRR4240389_cleaner.fastq (были удалены 12 688 056 первых строчек):

tail -n 12 688 056 SRR4240389_cleaner.fastq>half_cleaner.fastq

Далее была произведена подготовка k-меров (k=2) и их сборка по алгоритму, описанному выше с помощью следующих команд:

velveth half_velveth 29 -fastq -short half_cleaner2.fastq
velvetg half_velveth/

Ниже, в Таблице 5 приведены для сравнения характеристики двух проведенных сборок с k=29. Как видно, качетсво сборки существенно снизилось.

Таблица 5. Сравнение характиристик сборки всех ридов и половины ридов.

Число ридов	N50	Максимальная длина контига		Размер сборки
3172014	1194	7456	1106	610586
6344028	2960	15093	1150	654172

Как видно из Таблицы 5. качество сборки сущестенно снизилось, параметр N50 уменьшился более чем в два раза, в два раза снизилась длина максимального контига, число контигов при этом фактичесски не изменилось, мало изменился и размер собранного генома.

Главнaя страница