Этап 1. Геномная сборка

Организм: Rothia nasimurium

TaxID:85336 (взято из аннотации протеома в UniProt)

Код доступа GenBank: GCA_002087015.1

Источник: взят из записи протеома UP000192359 (поле «Genome assembly»).

Код доступа RefSeq: GCF_002087015.1

Источник: найден при поиске в NCBI Assembly (Datasets Genome) по коду GCA_002087015.1. Данная сборка имеет версию в RefSeq, которая и используется в дальнейшей работе.

Источник: страница протеома в базе UniProt Proteomes.

Этап 2. Скачивание последовательности генома и таблицы локальных особенностей (GFF3)

Команда для загрузки генома и GFF3:

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

Команда для распаковки архива:

unzip ncbi_dataset.zip

Этап 3. Поиск и трансляция открытых рамок считывания (ORF)

1. Определение генетического кода

Находим TaxID по сборке RefSeq:

esearch -db assembly -query GCF_002087015.1 | esummary | xtract -pattern DocumentSummary -element Taxid

Вывод: 85336

Проверяем номер генетического кода:

efetch -db taxonomy -id 85336 -format xml | grep -A1 'GCId'

Результат: 11

2. Поиск ORF и трансляция (getorf)

Задаём путь к файлу генома (получен на этапе 2):

GENOME_FILE="ncbi_dataset/data/GCF_002087015.1/GCF_002087015.1_ASM208701v1_genomic.fna"

Запускаем getorf для поиска рамок между стоп-кодонами (-find 0), минимальная длина 150 нуклеотидов (50 аминокислот), таблица 11:

getorf -sequence $GENOME_FILE -outseq orf_translations.faa -table 11 -minsize 150 -find 0

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

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

4. Проверка длины трансляций (infoseq)

Убеждаемся, что нет последовательностей короче 50 аминокислот:

infoseq -sequence orf_translations.faa -only -length | sort -n | head

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

Получаем последовательности трёх прокариотических ДНК-метилтрансфераз (P0AED9, P0AEE8, P23941) и сохраняем в query.fasta:

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

Этап 5. Поиск по сходству последовательностей (blastp)

Запускаем blastp с запросом query.fasta (этап 4) против базы proteome, выводим в формате 7:

blastp -query query.fasta -db proteome -outfmt 7 > blast_results.tsv

Ссылка: blast_results.tsv

Лучшая находка: NZ_LXWF01000007.1_268, координаты в геноме: 2071–3396 (обратная цепь), E-value: 1.14e-17, bit score: 82.8

Извлекаем CDS из GFF3 для того же контига:

grep "NZ_LXWF01000007.1_268" ncbi_dataset/data/GCF_002087015.1/genomic.gff | grep -w "CDS" | cut -f4,5,7,9 > CDS.tsv

Добавляем строку найденной ORF, сортируем по начальной координате и смотрим соседей (±3 записи):

echo -e  '2071\t3396\t-\tNZ_LXWF01000007.1_268' | cat - CDS.tsv | sort -n -k1,1 | grep -C 3 'NZ_LXWF01000007.1_268'> neighbors.tsv

Ссылка на neighbors.tsv: neighbors.tsv

Анализ neighbours.tsv показывает, что найденная ORF практически полностью перекрывается с аннотированной CDS WP_083090892.1 [2068 - 3393]. Аннотация CDS указывает на ДНК-цитозин-5-метилтрансферазу (ген dcm), что полностью соответствует цели поиска (гомолог P0AED9).

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

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