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

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

Команда:

fastqc chr14.1.fastq

Выдаёт информацию о качестве прочтений

Число чтений

18189

Длина чтений

41-51

Комментарий

Не было необходимости в использовании программы trimmomatic, поскольку как видно чтения имеют хорошее качество. Большого разброса по качеству не было, да и раницы не примыкали к красной зоне.

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

Команды:

1_ hisat2-build chr14.fasta chr14

Функция - индексирование референсной последовательности

2_ hisat2 --no-softclip -x chr14 -U chr14.1.fastq -S vir.sam

Функция - построение выравнивания референсной последовательности и прочтений

3_ samtools view vir.sam -b -o vir.bam

Функция - перевод выравнивания в бинарный формат

4_ samtools sort vir.bam -T file.txt -o virs.bam

Функция - сортировка выравнивания чтений с референсом по координате в референсе начала чтения

5_ samtools index virs.bam

Функция - индексирование отсортированногo .bam файла

Стоит отметить, что из команды hisat2 была убрана опция --no-spliced-alignment, так как в случае транскриптов, прошедших процессинг, некоторые участки изначальной последовательности ДНК уже вырезаны (поэтому при картировании нужно разрешить программе выравнивать фрагменты транскриптов по отдельности)

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

После вызова программы Hisat2 на экран выводится информация о выравнивании, в том числе количество выравниваний различных прочтений. Для данной реплики 47 прочтений не выровнялись, 18142 прочтения выровнялось 1 раз, 0 прочтений выровнялось > 1 раза.

4.Подсчёт чтений

В данном задании была использована программа htseq-count. Далее приведены некоторые её параметры:

-f - формат входных данных. Может использоваться формат .sam или .bam; по умолчанию стоит формат .sam

-s - зависимость данных от направления цепи. Yes - считывание с прямой цепи, no - нет зависимости от напрвления цепи, reverse - считывание с обратной цепи. По умолчанию стоит yes

-i - какой использовался атрибут GFF. По умолчанию стоит gene_id

-m - Режим обработки чтения, перекрывающий несколько функций. Возможными значениями являются объединение, строгое пересечение и пересечение - непустое; по умолчанию: объединение

Команда:

htseq-count -f bam -s no virs.bam /P/y14/term3/block4/SNP/rnaseq_reads/gencode.v19.chr_patch_hapl_sacaff.annotation.gtf >> itog.txt

Выдача скрипта:

__no_feature 27

__ambiguous 0

__too_low_aQual 0

__not_aligned 47

__alignment_not_unique 0

Далее была применена команда grep -w -v 0 > genes.count, в результате которых отсеялись строки с числом ридов равным 0. В результате высветился только один ген c 18115 ридами:

ENSG00000080824.14 18115

Некоторые чтения легли в неаннотированные участки генома.

Информация о гене:

A12M1 adenovirus-12 chromosome modification site 1C

Испльзование различных режимов обработки перекрывания ридов (union, intersection-strict, intersection-nonempty) не привело к изменениям.