Как и в практикуме 11, перед началом работы из директории /P/y14/term3/block4/SNP/rnaseq_reads в папку /nfs/srv/databases/ngs/s.kozyulina были скопированы файлы chr17.fasta, chr17.1.fastq и chr17.2.fastq.
Команды: fastqc chr17.1.fastq , fastqc chr17.2.fastq.
Выдача: файлы chr17.2_fastqc.html
и chr17.1_fastqc.html.
Характеристики | chr17.1.fastq | chr17.2.fastq |
График "Per base quality" из программы FastQC | ||
Число чтений | 10407 | 8195 |
Длина чтений | 25-51 | 32-51 |
Вывод:
Так как чтения имеют хорошее качество (quality score превышает 30), программа
trimmomatic не применялась.
№ команды | Команда | Результат выполнения команды |
1 | export PATH=${PATH}:/home/students/y06/anastaisha_w/hisat2-2.0.5 | вызывает программу hisat2-2.0.5, лежащую в данной директории |
2 | hisat2-build chr17.fasta chr17 | |
3 | hisat2 -x chr17 -U chr17.1.fastq --no-softclip > align1.sam
hisat2 -x chr17 -U chr17.2.fastq --no-softclip > align2.sam |
Строит выранивание прочтения и референса, сохраняет результаты в отдельный файл align1.sam, align2.sam |
В отличие от прак.11, из 3ей команды был убран параметр --no-spliced-alignment, так как в прак.12 мы работаем с РНК, в которых уже вырезаны интроны, поэтому относительно последовательности генома в них есть разрывы.
После вызова команды 3 на экран выводится следующая информация о выравнивании.
В первой реплике из 10407 прочтений 329 не были выравнены вовсе, 10001 - 1 раз, 77 - больше 1 раза.
Во второй реплике из 8195 прочтений 292 не были выравнены вовсе, 7847 - 1 раз, 56 - больше 1 раза.
№ команды | Команда | Результат выполнения команды |
4 | samtools view align1.sam -bo align1.bam
samtools view align2.sam -bo align2.bam |
Переводит выравнивание в бинарный формат |
5 | samtools sort align1.bam -T file.txt -o alignsort1.bam
samtools sort align2.bam -T file.txt -o alignsort2.bam |
Cортирует выравнивание чтений с референсом по координате в референсе начала чтения |
6 | samtools index alignsort1.bam
samtools index alignsort2.bam |
Индексирует отсортированный .bam файл |
№ команды | Команда | Результат выполнения команды |
1 | /P/y14/term3/block4/SNP/bedtools2/bin/bedtools bamtobed -i alignsort1.bam > chr17_1.bed | Вызывает программу bedtools из указанной директории, переводит файл формата .bam в формат .bed |
2 | /P/y14/term3/block4/SNP/bedtools2/bin/bedtools intersect
-a /P/y14/term3/block4/SNP/rnaseq_reads/gencode.genes.bed -b chr17_1.bed -c > sim.bed
sort -k 6 -r sim.bed > simsorted.bed |
Cоздается файл, где со списком генов идут их покрытия - количества совпадений с chr17_1.bed
(отображение покрытий обеспечивается параметром -с)
Перемещает гены с ненулевым покрытием в верх файла simsorted.bed |
3 | /P/y14/term3/block4/SNP/bedtools2/bin/bedtools intersect -a /P/y14/term3/block4/SNP/rnaseq_reads/gencode.genes.bed -b chr17_1.bed -u > chr17_1u.bed | Cоздается файл, содержащий строки лишь с генами, имеющими ненулевое покрытие (за счёт параметра -u) |
4 | /P/y14/term3/block4/SNP/bedtools2/bin/bedtools intersect -wa -wb -a chr17_1.bed -b /P/y14/term3/block4/SNP/rnaseq_reads/gencode.genes.bed > chr17_1_wa_wb.bed | Cоздается файл, где любые перекрывающиеся участки выводятся одной строкой, содержащей названия этих участков в файле chr5_1.bed и в разметке соответственно (за счёт параметров -wa и -wb) |
По итогам анализа выдачи программы bedtools были получены следующие данные:
Ген | Координаты | Функция | Количество ридов | Направление |
KPNB1 | 45726842 45762871 | protein_coding. Кодирует кариоферин-β1 - транспортный фактор, участвующий в переносе белков из цитоплазмы в клеточное ядро | 10063 | - |
RP11-580I16.2 | 45726762 45726786 | antisense | 1 | - |
1. Получите из файла в выравниванием файл с чтениями в формате fastq.
Команда: /P/y14/term3/block4/SNP/bedtools2/bin/bedtools bamtofastq -i alignsort1.bam -fq reads.fq
2. Получите файл с нуклеотидной последовательностью (.fasta) для одного из покрытых Вашими чтениями генов.
Команда: /P/y14/term3/block4/SNP/bedtools2/bin/bedtools getfasta -fi chr17.fasta -bed get_seq.bed > get_seq.fasta
3. Разбейте свою хромосому на фрагменты по 1 млн нуклеотидов. Какова длина хромосомы в нуклеотидах? Сколько в результате получилось интервалов?
Команда: bedtools makewindows -g chr17.fasta.fai -w 1000000 > mill.txt
Длина хромосомы в нуклеотидах - 81195210, количество интервалов - 82.