Практикум 11

Задача: Найти и описать полиморфизмы у пациента

Часть I: подготовка чтений

0. Создание рабочей директории.


В задании 0 требовалось В директории /nfs/srv/databases/ngs/ создайте свою директорию и скопируйте в нее Ваши файлы с ридами (.fastq) и хромосомой (.fasta).

1. Анализ качества чтений.


Сделайте контроль качества Ваших чтений с помощью программы FastQC. Комментарий: программа FastQC установлена на kodomo, её можно вызвать командой "fastqc file.fastq", где file.fastq — имя файла с чтениями. Версию с графическим интерфейсом можно поставить на свой компьютер. В результате работы программы Вы получите архив (.zip), который содержит отчет о программе в виде html файла.

Ссылка на отчет fastqc до очистки

Ссылка на отчет fastqc после очистки

2. Очистка чтений.

Очистка чтений проводилась с помощью программы Trimmomatic. Отрежьте с конца каждого чтения нуклеотиды с качеством ниже 20, оставьте только чтения длиной не меньше 50 нуклеотидов. Так как программа Trimmomatic установлена на kodomo, то первым делом вызываем ее: java -jar /usr/share/java/trimmomatic.jar SE -phred33 infile.fastq outfile.fastq step где infile.fastq и outfile.fastq — входной и выходной файлы с чтениями, а step — выражение, указывающее, какую операцию производить. Итак, команды: TRAILING: Cut bases off the end of a read, if below a threshold quality MINLEN: Drop the read if it is below a specified length Сделайте анализ качества очищенных чтений с помощью FastQC; сравните с прежней выдачей FastQC. java -jar /usr/share/java/trimmomatic.jar SE -phred33 infile.fastq outfile.fastq TRAILING:20 java -jar /usr/share/java/trimmomatic.jar SE -phred33 infile.fastq outfile.fastq MINLEN:50

Рис.1. "Per base quality" до очистки

Рис.2. "Per base quality" после очистки

Часть II: картирование чтений

3. Картирование чтений.


Откартируйте очищенные чтения с помощью программы Hisat2. Этапы Сначала необходимо проиндексировать референсную последовательность; команда hisat2-build Команда:export PATH=${PATH}:/home/students/y06/anastaisha_w/hisat2-2.0.5 hisat2-build chr21.fasta chr21 В результате возникли файлы chr21.1.ht2,...,chr21.6.ht2 Затем построить выравнивание прочтений и референса в формате .sam. Запустите hisat2 с параметрами --no-spliced-alignment и --no-softclip Команда: hisat2 -x chr21 -U chr21afterafter.fastq --no-spliced-alignment --no-softclip chr21.sam

4. Анализ выравнивания

Переведите выравнивание чтений с референсом в бинарный формат .bam. Используйте пакет samtools, команда view: samtools view; Команда: samtools view -b chr21.sam > chr21.bam Отсортируйте выравнивание чтений с референсом (получившийся после картирования .bam файл) по координате в референсе начала чтения; команда samtools sort; Команда: samtools sort chr21.bam chr21sort В итоге имеем файл chr21sort.bam Проиндексируйте отсортированный .bam файл командой samtools index Команда: samtools index chr21sort.bam Появился новый файл chr21sort.bam.bai

Рис.3. Информация из Hisat2

Часть III: Анализ SNP

5. Поиск SNP и инделей.

Создайте файл с полиморфизмами в формате .bcf; команда samtools mpileup -uf. Опции и формат описаны в руководстве. Команда: samtools mpileup -uf chr21.fasta chr21sort bam -o chr21poly.bcf. Создайте файл со списком отличий между референсом и чтениями в формате .vcf. Используйте команду "bcftools call -cv" пакета bcftools. Команда: bcftools call chr21poly.bcf -cv > chr21.vcf Найдите и опишите в отчете три полиморфизма из .vcf файла. Для каждого приведите: кординату; тип полиморфизма: замена, вставка или делеция; что было в референсе, что найдено в чтениях; глубина покрытия данного места; качество чтений в данном месте.

Координата Тип полиморфизма Референс Чтение Глубина покрытия Качество чтения
16334963ЗаменаТ G75225.009
16336804ЗаменаТ С31221.999
16340289ЗаменаС Т30169.009

6. Аннотация SNP.


С помощью программы annovar проаннотируйте только полученные snp (индели не надо!). Используйте базы данных: refgene, dbsnp, 1000 genomes, GWAS, Clinvar.
Комментарий: программа установлена на kodomo: /nfs/srv/databases/annovar. Для работы с программой annovar из .vcf файла необходимо получить файл, с которым умеет работать эта программа. Сделать это можно с помощью скрипта convert2annovar.pl. См. руководство.
Для аннотации файла с snp с помощью предложенных баз данных используйте скрипт annotate_variation.pl.

Команда: export PATH=${PATH}:/nfs/srv/databases/annovar perl convert2annovar.pl -format vcf4 chr21.vcf > chr21.avinput

Для аннотации файла с snp с помощью предложенных баз данных используйте скрипт annotate_variation.pl. В руководстве можно найти всю необходимую информацию о работе с программой. Например, узнать, какие из Ваших snp имеют rs, можно с помощью команды:

annotate_variation.pl -filter -out outputfile -build hg19 -dbtype snp138 inputfile.human humandb/
где inputfile.human — входной файл, полученный после обработки .vcf с помощью convert2annovar.pl (расширение не имеет значения); outputfile — выходной файл; humandb/ — директория, в которой лежат базы данных (все необходимые базы данных уже есть на kodomo, пользоваться опцией -downdb не надо!); snp138 — база данных, с которой вы работаете. Базы данных в annovar часто обновляются, для корректного запуска программы всегда нужно знать, какая версия какой базы данных у Вас скачена. Для вас: refgene — refGene; dbsnp — snp138; 1000 genomes — 1000g2014oct; GWAS — gwasCatalog; Clinvar — clinvar_20150629. В Annovar существуют 3 типа аннотаций по базам данных, основанных на: генной разметке (gene-based annotation); разметке других регионов генома (region-based annotation); фильтрации (filter-based annotation).

Команды, с помощью которых можно проаннотировать полиморфизмы по необходимым базам данных:

1. refgene - gene-based annotation

Команда: perl annotate_variation.pl -out refgene -build hg19 chr21.avinput /nfs/srv/databases/annovar/humandb/

Рис.3. Команда "refgene"

Результаты:

refgene.log

refgene.variant_function

refgene.exonic_variant_function

В результате получается таблица с указанием участков с заменой аминокислот. Таким образом имеем 3 синонимичных SNV и 1 несинонимичную замену. Далее в выдаче программы указывается координаты, и что именно на что заменилось. Синонимичные замены: C -> T (так и остался глицин); G -> C (так и остался серин); A-> G (так и остался пролин). Несинонимичная: A -> G (серин -> глицин). (это то, что было в файле refgene.exonic_variant_function). В другом выдаваемом файле описаны все замены, и не только в экзонах. Еще указаны гены, в которых это произошло: NRIP1, UBASH3A, AGPAT3. Также указано качество и глубина, гомо-/гетерозиготность. Глубина не очень хорошая, ее значние совсем невелико, в основном <5.

2.dbsnp - filter-based annotation

Команда: perl annotate_variation.pl -filter -out dbsnp -build hg19 -dbtype snp138 chr21.avinput /nfs/srv/databases/annovar/humandb/

Рис.4. Команда "dbsnp"

Результаты:

dbsnp.log

dbsnp.hg19_snp138_dropped

dbsnp.hg19_snp138_filtered

Здесь в результате можно увидеть сколько snp имеют rs - 63 в итоге.

3.1000 genomes - filter-based annotation

Команда: perl annotate_variation.pl -filter -dbtype 1000g2014oct_all -out 1000g -buildver hg19 chr21.avinput /nfs/srv/databases/annovar/humandb/

Рис.5. Команда "1000g"

Результаты:

1000g.log

1000g.hg19_ALL.sites.2014_10_dropped

1000g.hg19_ALL.sites.2014_10_filtered

Программа по аннотации по 1000 genomes позволяет оценить встречаемость данного snp в среднем у людей. Самая большая частота - 0.81869. Самая маленькая - 0.0131789.

4.Gwas - region-based annotation

Команда: perl annotate_variation.pl -regionanno -build hg19 -out ex1 -dbtype gwasCatalog example/ex1.avinput humandb/

Рис.6. Команда "gwas"

Результаты:

gwas.log

gwas.hg19_gwasCatalog

Gwas снабжает сведениями про ассоциированность с болезнями. В данном случае snp связаны со следубщими болезнями: Cognitive performance, Type 1 diabetes, Phospholipid levels (plasma).

5.Clinvar - filter-based annotation

perl annotate_variation.pl example/ex1.avinput humandb/ -filter -dbtype clinvar_20140211 -buildver hg19 -out ex1

Рис.7. Команда "clinvar"

Результаты:

clinvar.log

clinvar.hg19_clinvar_20150629_dropped

clinvar.hg19_clinvar_20150629_filtered


(С табличкой возникли проблемы, но она будет здесь и очень скоро.)

© Cherkashina Anastasia 2017