Поиск ДНК-метилтрансфераз в геноме Halorarum salinum средствами EMBOSS, EDirect, NCBI Datasets CLI и blast+

Для проведения работы была выбрана архея Halorarum salinum. Белковые записи её протеома были получены с помощью поиска в UniProt Proteomes. Через запись UniProt Proteomes, соответствующей данному протеому, были получены коды доступа соответствующей ему геномной сборки.

Идентификатор протеома в UniProt Proteomes: UP000509626

Код доступа геномной сборки GenBank: GCA_013402875.1

Код доступа геномной сборки RefSeq: GCF_013402875.1

Страница сборки в Datasets Genome

Файлы, содержащие геномную последовательность и таблицу особенностей в формате gff3, были получены с помощью следующей команды:

datasets download genome accession GCF_013402875.1 --include genome,gff3

Скачанный архив был распакован с помощью команды unzip:

unzip ncbi_dataset.zip

Для поиска ORF необходимо было определить, какой генетический код характерен для H. salinum. Эта информация содержится в записи NCBI Тaxonomy, соответствующей археи. Данная запись была получена в формате xml с помощью следующей команды:

esearch -db 'taxonomy' -query '"halorarum salinum"[Scientific Name]' | efetch -format 'xml' >tax.txt

Информация об используемом генетическом коде содержится в поле "GCid". Это таблица №11.

Далее был проведён поиск ORF в геномной последовательности H. salinum с помощью команды getorf из пакета EMBOSS:

getorf GCF_013402875.1_ASM1340287v1_genomic.fna orfs.fasta -table 11 -minsize 150

Указанные параметры задают таблицу №11 генетического кода для чтения, минимальную длину (в нуклеотидах) выводимых рамок считывания, равную 150. Были получены рамки между стоп-кодонами (параметр по умолчанию) и записаны в файл 'orfs.fasta'. По полученным последовательностям была создана локальная база данных для blastp, названная proteome:

makeblastdb -in orfs.fasta -dbtype 'prot' -out proteome

Примечание: с помощью следующей команды выдача команды getorf, указанной ранее, была проверена на наличие результатов, имеющих длину менее 50 аминокислот (150 нуклеотидов):

infoseq -filter orfs.fasta -only -length | sort -n |& less

Параметры -only -length для команды infoseq в сочетании позволяют вывести лишь длины последовательностей, что удобно для проверки.

По результатам проверки таких последовательностей обнаружено не было, а значит команда была задана корректно.

Для поиска ДНК-метилтрансфераз в геноме H. salinum по сходству последовательностей с другими ДНК-метилтрансферазами, были получены последовательности ДНК-метилтрансфераз Bacillus amyloliquefaciens (Swiss-Prot AC: P23941) и E. coli (P0AED9, P0AEE8). Для этого использовалась команда:

echo 'sw:P0AED9 sw:P0AEE8 sw:P23941' | tr " " "\n" | seqret -filter @stdin >query.fasta

Последовательности были сохранены в файл 'query.fasta'.

С помощью локального blastp было проведено выравнивание данных трёх последовательностей с транслированными ORF из генома H. salinum (то есть выравнивание проводилось по созданной ранее локальной базе proteome). Применённая для этого команда:

blastp -query query.fasta -db proteome -outfmt 7 -out res.tab

Табличная выдача blastp

Лучшая по весу (52.4 bits) находка - NZ_CP058579.1_11181. Её координаты в геноме - 2219157-2220395. Судя по результатам выравнивания, последовательность данной находки гомологична m6A-метилазе E. coli (Swiss-Prot AC: P0AEE8): в данном выравнивании есть достаточно достаточно продолжительные идентичные и схожие участки.

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

<genomic.gff cut -f 3,4,5,7,9 | grep "^CDS" | cut -f 2,3,4,5 >CDS.tsv

С помощью этого файла был проведён поиск CDS, расположенных рядом с найденной ORF. К строкам из CDS.tsv была добавлена строка, содержащая информацию об искомой ORF. Затем все строки были отсортированы по первому столбцу - меньшей из координат последовательностей. Десять строк, оказавшихся рядом с найденной ORF, были записаны в файл 'neighbors.tsv'. Это было реализовано с помощью команды:

echo -e "2219157\t2220395\t+\tFOUND_ORF" | cat - CDS.tsv | sort -k 1 -n | grep -C 5 "FOUND_ORF" > neighbors.tsv

По содержимому 'neighbors.tsv' было определено, что с рассматриваемой ORF пересекаются 2 CDS, одна из которых полностью лежит в данной ORF (вторая пересекается с ORF на незначительном участке и задана в другой рамке считывания). Координаты этой CDS - 2219346-2220398. Это означает, что её стоп-кодон является стоп-кодоном, ограничивающим рассматриваемую ORF. В аннотациях, приведенных в таблице особенностей, видим, что данная CDS соответствует аденин-метилтрансферазе (product=DNA adenine methylase), а значит, нам действительно удалось найти ДНК метилтрансферазу в геноме H. salinum (а именно m6A-МТазу). Код доступа данной CDS - WP_246308032.1.

Проверим, возможно ли было найти данную CDS прямым поиском по аннотациям. Для этого использовалась следующая команда:

elink -db 'nuccore' -id 'NZ_CP058579.1' -target 'protein' | efilter -query '2.1.1.72[EC/RN Number]' | efetch -format 'acc'

Здесь NZ_CP058579.1 - код доступа нуклеотидной последовательнности из базы nucleotide, в которой была обнаружена интересующая нас ORF. Т.к. поиск по названию фермента в NCBI проводить сложно, для поиска использовался EC код аденин-метилтрансферазы - 2.1.1.72.

В результате проведения поиска была обнаружена лишь одна запись с кодом доступа WP_246308032.1. Эта та же запись, которую удалось найти через транслирование ORF геномной последовательности H. salinum.