EMBOSS, Entrez Direct, NCBI Datasets, BLAST+
В данном практикуме предлагалось с помощью инструментов EMBOSS, EDirect, NCBI Datasets CLI и BLAST+ найти одну из ДНК-метилтрансфераз выбранного в первом семестре прокариота по последовательности и по аннотации.
Получение АС геномной сборки и TaxID организма
Исследуемый прокариот — архея Natronomonas pharaonis DSM 2160. По идентификатору протеома в UniProt Proteomes была найдена соответствующая запись, откуда была извлечена информация о коде доступа (Accession code, AC) сборки генома в GenBank и RefSeq. Со страницы сборки в NCBI Genome была получена информация о версии сборки в RefSeq. Все необходимые идентификаторы, коды доступа и ссылки приведены ниже.
1. Идентификатор протеома: UP000002698
2. Ссылка на запись из UniProt Proteomes: https://www.uniprot.org/proteomes/UP000002698
3. Код доступа сборки GenBank: GCA_000026045.1
4. Ссылка на страницу сборки в NCBI Genome: https://www.ncbi.nlm.nih.gov/datasets/genome/GCF_000026045.1/
5. Код доступа сборки RefSeq: GCF_000026045.1
Далее с помощью приведённых ниже команд был скачен и распокован архив, содержащий файлы с полной последовательностью генома
и с таблицей его локальных особенностей в формате gff3.
unzip ncbi_dataset.zip
Поиск и трансляция открытых рамок считывания
Рассматриваемый организм Natronomonas pharaonis DSM 2160 использует таблицу генетического кода №11, что было выявлено после получения записи про таксон из базы NCBI Taxonomy (тег GeneticCode, внутри него тег GCId, внутри 11):
где -id — идентификатор организма в базе NCBI Taxonomy (Taxonomy ID), который был получен из поля Taxon ID в записи протеома в базе UniProt Proteomes.
Затем были получены трансляции открытых рамок считывания с помощью программы getorf пакета программ EMBOSS. Рамки выбирались между стоп-кодонами, чтобы не потерять фрагмент последовательности из-за неверного старт-кодона, при этом их трансляции должны были быть не короче 50 аминокислотных остатков, с целью отбросить маловероятные короткие фрагменты. Результат был сохранён в промежуточном файле proteome.orf. Затем полученный файл использовался в качестве белковой базы для blastp (файл proteome), которая создавалась с помощью соответствующей команды makeblastdb. В конце с помощью программы infoseq было проверено, что среди полученных трансляций нет тех, которые короче 50 а/о.
getorf -sequence ncbi_dataset/data/GCA_000026045.1/GCA_000026045.1_ASM2604v1_genomic.fna -table 11 -minsize 150 -find 0 -outseq proteome.orf
makeblastdb -in proteome.orf -dbtype prot -out proteome
infoseq -only -length proteome.orf | sort -n | less
Получение последовательностей гомологичных метилтрансфераз
Изначально в файл term3/pr9/tmp_MT.txt были записаны USA соответствующих метилтрансфераз: sw:P0AED9 (m5C), sw:P0AEE8 (m6A), sw:P23941 (m4C). Далее с помощью команды seqret были получены все три белковые последовательности из базы данных Swiss-Prot. Результат был записан в файл query.fasta.
Поиск по сходству последовательностей
С помощью локального blastp по белковой базе proteome, используя содержимое query.fasta в качестве запроса, проводился поиск по сходству последовательностей гомологичных метилтрансфераз. Затем рассматривалась только лучшая по весу находка — с идентификатором CR936257.1_23884, гомологичная m4C-метилтрансферазе, координатами в геноме: 1144506 - 1142851 (были получены из файла proteome.orf при поиске по идентификатору находки), весом 59.3 бита. Сомнений по поводу гомологичности находки не возникло, так как evalue = 6.6*10^-10 (что говорит о высокой достоверности находки). Ссылка на выдачу.
Далее с помощью соответствующих конвейеров на Bash было выяснено, какие CDS из таблицы локальных особенностей генома рассматриваемого организма располагаются рядом с находкой. Всего таких последовательностей оказалось 7 штук (2 на "+"-цепи и 5 на "-"-цепи), одна из которых пересекается с найденно открытой рамкой считывания и кодирует продукт site-specific DNA-methyltransferase. Координаты почти совпадают, следовательно эта CDS и является гомологом m4C-метилтрансферазы. Ссылка на файл с выдачей. Конвейеры в Bash:
blastp -db proteome -query query.fasta -outfmt 7 > tmp_query_MT_fmt7.txt
grep "CR936257.1" ncbi_dataset/data/GCA_000026045.1/genomic.gff | grep "CDS" | cut -f4,5,7,9 > CDS.tsv
echo -e '1144506\1142851\t-\tFOUND-ORF' | cat - CDS.tsv | sort -n | grep -C 3 'FOUND-ORF' > neighbors.tsv
Поиск по аннотациям кодирующих участков
В завершение проводился поиск CDS по аннотациям кодирующих участков в геноме. Конвейер EDirect: