Практикум 9: Поиск ДНК-метилтрансфераз в геноме прокариоты

Глобальная задача: найти одну из ДНК-метилтрансфераз (m5C, m6A или m4C) в геноме бактерии/археи с использованием средств EMBOSS, EDirect, NCBI Datasets CLI и blast+.

Этап 0: Протеом бактерии

Использован протеом организма из практикума 8 семестра.

Команда для скачивания:

wget "https://www.uniprot.org/uniprot/?query=proteome:UP000001234&format=swiss" -O uniprot-proteome.swiss

Этап 1: Получение AC геномной сборки и TaxID организма

поиск TaxID:

поиск происходил с помощью конвееров на языке bash в несколько этапов:

Этап 2: Поиск сборки

datasets summary genome taxon 1616823 --as-json-lines | dataformat tsv genome

Команда выдала стену текста, сложную для осмысления. Для упрощения добавила в конвеер

cut -f1,2,3 | sort -u"

Это показало, что сборка одна и у неё есть 2 версии:

Для дальнейшей работы я выбрала версию RefSeq (так как она лучше аннотированна), AC сборки: GCF_000948415.1

Этап 3: Скачивание последовательности генома и таблицы локальных особенностей

Имея AC сборки, можно скачать нуклеотидную последовательность генома и таблицу локальных особенностей

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

Скачанный архив был в формате .zip

unzip ncbi_dataset.zip

Этап 4: Поиск и трансляция открытых рамок считывания

Прежде, чем искать открытые рамки считывания, нужно определить, какой вариант генетического кода использует организм. Эту информацию можно получить из записи про таксон в базе NCBI Тaxonomy, скачав её с помощью efetch в формате xml.

efetch -db taxonomy -id 1616823 -format xml

Генетический код: 11

Поиск ORF с помощью getorf

getorf -sequence ncbi_dataset/data/GCF_000948415.1/*.fna -outseq orf_translations.fasta -find 0 -minsize 150 -table 11

Пояснение:

Проверка:

infoseq orf_translations.fasta

Создание BLAST-базы:

makeblastdb -in orf_translations.fasta -dbtype prot -out proteome

Пояснение:

Программа созадала несколько файлов:

Этап 5: Получение последовательностей гомологичных метилтрансфераз

Все известные ДНК-метилтрансферазы прокариот, по всей видимости, содержат гомологичные каталитические домены. Однако они могут быть насколько далеки друг от друга, что сходство последовательностей стандартными средствами обнаружить не удается, но мы попробуем.

echo -e 'sw:P0AED9\nsw:P0AEE8\nsw:P23941' > list.txt

Надо скачать последовательности гомологичных метилтрансфераз (по сходству с которыми можно было бы осуществить поиск генов этих белков у данного организма), а именно:

Наиболее коротким способом, это можно сделать в 2 этапа: создать файл-список со USA, а далее скачать все последовательности из файла-списока

echo -e 'sw:P0AED9\nsw:P0AEE8\nsw:P23941' > list.txt

seqret @list.txt -outseq query.fasta

Этап 6: Поиск по сходству последовательностей

Наконец-то производим поиск гена метилтрансферазы у данного организма с помощью blastp, созданной ранее базы данных ORFs и скачанных в прошлом пункте последовательностей в качестве запроса:

Команда BLASTp:

blastp -query query.fasta -db proteome -out blast_results.txt -outfmt 7 

Результаты BLAST: blast_results.txt

Анализ лучшего хита:

Самая лучшая находка по весу, как можно видеть из выдачи, имеет идентификатор NZ_JYII01000006.1_3165

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

grep "NZ_JYII01000006.1" genomic.gff | grep "CDS" | cut -f4,5,7,9 > CDS.tsv

Файл срздали, теперь ищем сходства.

 echo -e '50\t295\t+\tFOUND-ORF' | cat - CDS.tsv | sort -n | grep -C 3 'FOUND-ORF' > neighbors.tsv

Файл с выводом: neighbors.tsv

Вывод: у ORF нет пересечения с CDS

Этап 6: Поиск по аннотациям кодирующих участков

Поскольку лучшая находка - m4C-метилтрансфераза, будем искать по EC-коду 2.1.1.113.

esearch -db nuccore -query "GCF_000948415.1" | elink -target protein | efilter -query "(2.1.1.37 OR 2.1.1.72 OR 2.1.1.113)" | efetch -format 'fasta'

В результате найдены две аннотированные ДНК-метилтрансферазы: DNA adenine methylase (WP_235285025.1, m6A-тип) и DNA cytosine methyltransferase (WP_044828408.1, m5C-тип). Однако m4C-метилтрансферазы в аннотации генома не обнаружено.

BLAST нашел ген метилтрансферазы, который отсутствует в официальной аннотации генома. Поиск по аннотациям обнаружил только два гена метилтрансфераз - m5C и m6A типа, а наш ген m4C типа в списке не значился. Это показывает, что BLAST является более мощным инструментом поиска, чем простой просмотр аннотаций. BLAST может находить гены по сходству последовательностей, даже если они не указаны в описании генома.