Одноконцевые RNA-Seq чтения из биологической реплики 1 (файл с чтениями chr8.1.fastq) картировались на референсную последовательность хромосомы 8 H. sapiens.
На основании разметки человеческого генома Gencode 19 для сборки hg19 определялись гены, с которых были транскрибированы отсеквенированные РНК выбранной реплики.

Все полученные файлы хранятся в рабочей директории /nfs/srv/databases/ngs/stepan_puhov/pr12 в соответствующих поддиректориях.

Далее представлена таблица, отражающая последовательность действий и использованных команд (для точности описания указывается, в каких директориях производился запуск соответствующих команд):

С помощью программы FastQC был получен отчёт о качестве чтений.

Из раздела отчёта "Basic statistics" видно, что в исследуемом образце всего 17763 чтений, длина чтений составляет от 25 nt до 51 nt, ни одно чтение не помечено как "чтение плохого качества" и используется кодировка качества Phred+33.
Из раздела отчёта "Sequence length distribution" понятно, что большинство чтений имеют длину более 49 nt:

Картирование чтений на референсную последовательность было произведено программой hisat2 c отключенной опцией --no-spliced-alignment.
Отключение этой опции позволяет программе разделять чтение на участки и картировать полученные участки независимо на одну и ту же цепь референса в предположении о наличии "сплайсинг-сайта" (вырезанного интрона) в данном чтении. Для чтения, картированного таким способом, в поле CIGAR sam-файла указывается длина прерывающего его интрона (кодом в регулярном выражении [0-9]+N); с помощью программы grep было подсчитанно, что так картировалось 4514 чтений.

Из отчёта картирования программы hisat2 видно, что из 17763 данных программе чтений 304 чтения не картировались, а из остальных чтений 17455 картировались ровно 1 раз и 4 чтения картировались более 1 раза. Таким образом, всего картировалось 98.29% чтений.
С помощью программы grep было выяснено, что у всех 17455 чтений, что картировались лишь 1 раз, было высокое качество картирования (Q = 60), а у 4 чтений, которые картировались более 1 раза, качество картирования было очень низким (Q = 1 или Q = 0). Можно заключить, что картирование, в целом, имеет высокое качество.

Для дальнейшей работы файл sam-формата с картированными чтениями с помощью программ пакета Samtools был переведён в bam-формат, а полученный bam-файл был отсортирован по координате картирования чтений и проиндексирован.

Подсчёт чтений, картированных на гены, определяемые согласно разметке Gencode 19 сборки hg19, производился программой htseq-count.

Программа принимает файл с картированными чтениями и файл формата gtf или gff с разметкой генома и для каждого участка генома, выделенного данной разметкой, (так называемого feature) считает число попавших на него чтений по 1 из 3 алгоритмов: union (для каждого чтения учитывается вхождение в участок генома (feature) хотя бы 1 нуклеотидом), intersection-strict (для каждого чтения учитывается вхождения в участок генома (feature) только по всем нуклеотидам чтения) или intersection-nonempty (для каждого чтения учитывется вхождение в feature по всем нуклеотидам чтения, входящих хотя бы в какой-нибудь feature).
В данном случае за features принимались гены, как совокупность их экзонов (то есть вхождение чтения в feature по данному нуклеотиду − это вхождение этого нуклеотида чтения в какой-нибудь экзон данного гена). При этом 3 алгоритма подсчёта чтений могут давать следующие результаты в разных ситуациях:

• -f определяет формат входного файла с картированными чтениями
• -s значение "yes" считает чтения для feature, только если они на той же цепи, что и feature, значение "reverse" − только если чтения на противоположной цепи, значение "no" позволяет считать чтения, картирующиеся на обе цепи по отношению к данному feature
• -t определяет название для выбранного типа feature во входном gtf файле (default: exon)
• -i определяет аттрибут в записях gtf файла, значение которого принимается за feature ID (default: gene_id)
• -m определяет 1 из 3 способов подсчёта чтений (default: union)

ENSG00000104738.12      294
ENSG00000253729.3       15941
__no_feature    	1219
__ambiguous     	1
__too_low_aQual 	0
__not_aligned   	304
__alignment_not_unique  8

ENSG00000104738.12    	273
ENSG00000253729.3     	15511
__no_feature  		1671
__ambiguous   		0
__too_low_aQual       	0
__not_aligned 		304
__alignment_not_unique  8

• 304 чтения не были картированы (они были определены ещё программой hisat2);
• 8 выравниваний не были уникальными (это выравнивания для тех самых 4 чтений с более, чем 1 картированием; у каждого из них есть 2 координаты картирования, что видно из sam файла);
• 1 чтение частично лежит на пересечении экзонов 2 разных генов;
• 1671 чтение хотя бы частично попадают на интроны, а частично могут попадать на экзоны гена PRKDC или гена MCM4, из них 1219 чтений полностью попадают на интроны;
• 273 и 15511 чтений полностью попадают на экзоны генов MCM4 и PRKDC соответсвенно.

• PRKDC − ген каталитической суб'единицы ДНК-активируемой протеинкиназы, серин/треониновой протеинкиназы, участвующей в NHEJ-репарации.
• MCM4 − один из 8 паралогичных генов, белки которых образуют гетерооктамерный комплекс, выполняющий роль главной хеликазы при репликации ядерной ДНК.