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

Идентификатор протеома: UP000236647

Ссылка на страницу Proteomes

Kод доступа геномной сборки из GenBank GCA_002900385.1, версия прокариотической сборки в RefSeq GCF_002900385.1

Cсылкa на страницу сборки в Genome

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

Команда, с помощью которой была загружена последовательность и таблица локальных особенностей:

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

Команда, с помощью которой мы распаковываем архив с файлами:

unzip ncbi_dataset.zip

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

Находим TaxID по Accession:

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

Получили 1348

Проверяем, использует ли наша бактерию Таблицу 11 с помощью команды:

efetch -db taxonomy -id 1348 -format xml | grep -A 1 'GCId'

Получили 11, а это значит, что наша бактерия использует Таблицу 11

Задаем переменную с путем к нужному файлу:

GENOME_FILE="ncbi_dataset/data/GCF_002900385.1/GCF_002900385.1_ASM290038v1_genomic.fna"

Проверяем, нет ли трансляции короче, чем 50 аминокислот:

getorf -sequence $GENOME_FILE -outseq proteome.fasta -table 11 -minlen 50

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

Получаем последовательности и сохраняем их в query.fasta:

efetch -db protein -id P0AED9,P0AEE8,P23941 -format fasta > query.fasta

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

Создаем базу данных BLAST

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

Поиск в новой базе:

blastp -query query.fasta -db proteome_db -outfmt 7 -max_target_seqs 1 > blast_results.tsv

Лучшая находка и по весу, и по e-value это гомолог m5C-MTазы (P0AED9), который соответствует рамке считывания NZ_CP025420.1_40840. Координаты рамки [353-448]

Ищем строку заголовка в файле proteome.fasta, чтобы узнать полные координаты:

grep "NZ_CP025420.1_40840" proteome.fasta

Координаты рамки [1888843 - 1887032]

GFF_FILE="ncbi_dataset/data/GCF_002900385.1/genomic.gff"

Фильтруем GFF3 по контигу NZ_CP025420.1 и типу CDS, извлекаем столбцы 4, 5, 7, 9:

grep "^NZ_CP025420.1" $GFF_FILE | grep -w "CDS" | cut -f 4,5,7,9 > CDS.tsv

Используем код для поиска соседей, подставляя наши координаты:

echo -e '1888843\t1887032\t+\tNZ_CP025420.1_40840' | cat - CDS.tsv | sort -n -k1,1 | grep -C 3 'NZ_CP025420.1_40840' > neighbors.tsv
Скачать neighbors.tsv

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

esearch -db nuccore -query "GCF_002900385.1" | elink -target protein | efilter -query "2.1.1.37[ECNO]" | wc -l

Найдено 7 белков