Однонуклеотидные полиморфизмы, индели и сборка

Семестры • Третий семестр • Однонуклеотидные полиморфизмы, индели и сборка

Поиск однонуклеотидных полиморфизмов и инделей

В этом задании мне нужно было с помощью программ samtools и bcftools получить список однонуклеотидных полиморфизмов (SNP) и инделей (делеций и инсерций) для чтений, картированных на геномы хлоропласта и митохондрии из предыдущего задания. Программа bcftools принимает на вход файл в формате .bcf. Использованные команды:

#Создание файлов в формате .bcf
#-g: рассчитывает сходство генотипов, выдает файл в формате .bcf
#-u: выдает файл в формате .bcf (эта опция отчасти копирует первую)
#-f: принимает дополнительно файл с референсной последовательностью
samtools mpileup -ugf chloroplast.fasta chlor-mem.sorted.bam > chlor.bcf
samtools mpileup -ugf mitochondrion.fasta mito-mem.sorted.bam > mito.bcf

#Анализ файлов и построение списка полиморфизмов и инделей
bcftools view -vcg chlor.bcf > chlor.vcf
bcftools view -vcg mito.bcf > mito.vcf  

#Подсчет отдельно инделей и полиморфизмов                                                                         
grep 'INDEL;' chlor.vcf | wc -l                                         
grep 'DP=' chlor.vcf | wc -l                                            
grep 'INDEL;' mito.vcf | wc -l                                          
grep 'DP=' mito.vcf | wc -l                                             

Результаты поиска представлены в таблице 1. Файлы можно скачать: chlor.vcf, mito.vcf.

Сборка хлоропласта и митохондрии

Сборка геномов из очищенных ридов осуществляется программой velvet с использованием графов Де Брейна. Для сборки сначала надо создать банк k-меров программой velveth, затем собрать риды в контиги программой velvetg. От длины k-мера зависит параметр n50, показывающий минимальную длину контигов, покрывающих 50% последовательности. Максимальная поддерживаемая длина k-мера 31, длины должны быть нечетными. Я провел сборку для длин 31, 29, 27, 25, 23, 21. Наилучший результат оказался при длине 23:

Final graph has 428300 hodes and n50 of 238, max 5086, total 47119351, using 0/3892273 reads

То есть n50 составило 238, 50% последовательности покрываются контигами длины 238 и более.

10 самых длинных контигов описаны в таблице 2. Файл с ними можно скачать.

Чтобы понять, какой контиг соответствует какой органелле, я провел локальный бласт против генома митохондрии и генома хлоропласта. Первые девять контигов отлично (e-value 0.0) выравниваются с геномом митохондрии, но один контиг (ID 11229) не соответствует ни геному митохондрии, ни геному хлоропласта (с геномом хлоропласта не выравнился ни один контиг из топ-10 по длине). Скорее всего этот контиг собран неверно.

Файл с выдачей BLAST против митохондрии можно скачать.