1.Получение AC геномной сборки и TaxID организма
Как и в прошлом семестре, буду работать с бактерией Brucella canis ATCC 23365. Для референсного протеома UP000001385 код доступа в GenBank - GCF_000018525.1, в RefSeq - GCF_000018525.1. В базе данных NCBI Genomes найдена соответствующая ссылка на сборку.
TaxID исследуемой бактерии - 483179.
2. Скачивание последовательности генома и таблицы локальных особенностей
Чтобы скачать последовательности генома и таблицы локальных особенностей сборки GCF_000018525.1, воспользуемся командой:
datasets download genome accession GCF_000018525.1 --include gff3,genome
Пояснения:
--include gff3,genome — получить файл с последовательностью генома и таблицу особенностей (feature table) в формате gff3
Архив был разархивирован с помощью следующей команды:
unzip ncbi_dataset.zip
3. Поиск и трансляция открытых рамок считывания
Для того, чтобы проводить поиск открытых рамок считывания (ORF), нужно узнать, какой вариант генетического кода используется исследуемой бактерией. С помощью этой команды была скачана запись из базы данных NCBI Taxonomy:
efetch -db 'taxonomy' -id '483179' -format 'xml'
Пояснения:
-db — база данных, в которой проводится поиск
-id — TaxID бактерии
-format — формат выходного файла (xml)
Информация об используемой бактерией генетической таблице находится в поле GeneticCode. Она использует 11 таблицу - самый распространенный вариант генетического кода для бактерий, архей и пластид.
Команда для получения последовательностей открытых рамок считывания:
getorf -sequence ncbi_dataset/data/GCF_000018525.1/GCF_000018525.1_ASM1700v1_genomic.fna -outseq ORFs.fasta -table 11 -minsize 150
Пояснения:
-sequence — название файла генома бактерии
-outseq — название выходного файла, в котором будут записаны все ORFs
-table — номер таблицы генетического кода, используемой бактерией (выяснили, что 11)
-minsize 150 — в ORFs должно быть не менее 50 аминокислотных остатков (150 нуклеотидов = 50 АК)
Для того, чтобы проверить, что во всех ORFs находятся не менее 50 аминокислотных остатков, был составлен следующий конвейер:
infoseq ORFs.fasta -only -length | sort -n | head
Пояснения:
-only -length — нужно вывести только длину последовательностей
sort -n сортирует длины по убыванию, head выводит только первые несколько строчек. Минимальная длина оказалось равной 50, значит конвейер работает правильно.
Команда для создания базы белков для blastp:
makeblastdb -in ORFs.fasta -dbtype prot -out proteome
Пояснения
-in ORFs.fasta — название файла со всеми ORF длиной не менее 50 АК
-dbtype prot — белковая база данных
-out proteome — имя базы данных
4. Получение последовательностей гомологичных метилтрансфераз
Записать в файл коды доступа swissprot метилтрансфераз (P0AED9, P0AEE8, P23941) можно при помощи одного конвейера:
echo -e 'sw:P0AED9\nsw:P0AEE8\nsw:P23941' | seqret -filter @stdin query.fasta
С помощью ранее составленной базы данных proteome и последовательностей query.fasta была найдена ORF, гомологичная метилтрансферазе.
blastp -query query.fasta -db proteome -out methyltr_blastp.out -outfmt 7
Пояснения:
-query query.fasta — имя файла с последовательностями метилтрансфераз
-db proteome — ранее созданная мной база белков из ORF
-out methyltr_blastp.out — имя выходного файла
-outfmt 7 — табличная выдача в формате 7
Лучшая находка с весом выравнивания (91.7) и e-value = 1.31e-20 имеет идентификатор NC_010103.1_2675. Она гомологична P23941 (m4C-МТаза, Bacillus amyloliquefaciens).
Чтобы найти координаты NC_010103.1_2675 в геноме бактерии, составлен конвейер:
cat ORFs.fasta | grep 'NC_010103.1_2675'
Координаты: [490390 - 491622].
Далее был проведён поиск близких последовательностей по координатам. Для этого нужно было составить таблицу, которая содержит CDS той же последовательности (NC_010103.1), что и находка. Был использован конвейер:
cat ncbi_dataset/data/GCF_000018525.1/genomic.gff | grep 'NC_010103.1' | grep 'CDS' | cut -f4,5,7,9 > CDS_NC_010103.1.tsv
Затем был создан файл neighbors.tsv, в котором находятся близкие по координатам последовательности:
echo -e '490390\t491622\t-\tFOUND-ORF' | cat - CDS_NC_010103_1.tsv | sort -n | grep -C 3 'FOUND-ORF' > neighbors.tsv
Ссылка на файл neighbors.tsv
Из полученной таблицы видно, что есть последовательность, координаты которой практически полностью совпадают с находкой (490492-491625). Она принадлежит записи WP_012090866.1. В описании данной записи указано product=site-specific DNA-methyltransferase, поэтому, скорее всего, эта последовательность действительно является гомологом.
6. Поиск по аннотациям кодирующих участков
Далее нужно было выяснить, нашлась ли CDS, соответствующая моей находке, с помощью поиска по аннотации кодирующих участков в геноме. Для этого был составлен следующий конвейер:
elink -db nuccore -target protein -id 'NC_010103.1' | efilter -query 'EC[ECNO]' | esummary
Пояснения:
elink -db nuccore -target protein -id 'NC_010103.1' — ищет в базе данных Nucleotides запись с идентификатором NC_010103.1 и находит все соответствующие белковые последовательности в базе данных Protein.
efilter -query 'EC'[ECNO] — поиск соответствующих метилтрансфераз
esummary — 'синоним' efetch c опцией -format 'docsum', для более удобного просмотра полученных данных
Вместо EC я указывала 2.1.1.37 (у m5C), 2.1.1.72 (у m6A) и 2.1.1.113 (у m4C). Но не нашлось ни одной находки, соответствующей этим метилтрансферазам.
С помощью скрипта захотелось проверить, какие метилтрансферазы (EC - 2.1.1) закодированы на 1 хромосоме Brucella canis:
elink -db nuccore -target protein -id 'NC_010103.1' | efilter -query '2.1.1[ECNO]' | esummary | xtract -pattern DocumentSummary -element AccessionVersion Title
