Ранее я сравнивала протеом Streptomyces vietnamensis с протеомом Streptomyces coelicolor в качестве контроля, соотвественно требовалось исследовать на ДНК-метилтрансферазы геном S. vietnamensis.
Для этого сначала был найден TaxID этой бактерии в swiss файле с белками протеома бактерии с помощью данной команды:
zcat UP000031774.swiss.gz | grep '^OX' | cut -b17-23 | sort | uniq -c
Выдача состояла из одной строки с символами 362257.
Далее был получен список сборок, соответствующих таксону:
datasets summary genome taxon 362257 --as-json-lines | dataformat tsv genome | vd
Я выбрала геномную сборку с AC GCF_000830005.1, которая является RefSeq версией и референсной.
По AC геномной сборки выше скачала последовательность генома (GCF_000830005.1_ASM83000v1_genomic.fna) и таблицу локальных особенностей (genomic.gff), использовав команду ниже:
datasets download genome accession GCF_000830005.1 --include 'genome','gff3'
Разархивировала ncbi_dataset с искомыми файлами так:
unzip ncbi_dataset.zip
Перед поиском открытых рамок считывания необходимо было проверить, какая таблица генетического кода соответствует S. vietnamensis. Для этого я скачала запись про таксон с id 362257 из NCBI Тaxonomy:
efetch -db taxonomy -id "362257" -format xml > tax_s_vietnamensis.xml
Там были найдены поля "GeneticCode" и "GCId", где была указана таблица - 11.
После данного подтверждения можно было приступать к нахождению открытых рамок считывания:
getorf -sequence "ncbi_dataset/data/GCF_000830005.1/GCF_000830005.1_ASM83000v1_genomic.fna" -outseq orfs_s_vietnamensis.fasta -minsize 150 -table 11 -find 0
Конвейером с программой infoseq я проверила, что все последовательности из полученного файла orfs_s_vietnamensis.fasta состоят хотя бы из 50 а.о.:
infoseq -only -length -sequence orfs_s_vietnamensis.fasta | sort -n -u | head -n4
По тем же последовательностям из файла orfs_s_vietnamensis.fasta создала белковую базу, назвав её ORFs:
makeblastdb -in orfs_s_vietnamensis.fasta -dbtype prot -out ORFs
Найти последовательность ДНК-метилтрансферазы в геноме было предложено по сходству с последовательностями ДНК (цитозин-5)-метилтрансферазы E. coli (код доступа Swiss-Prot: P0AED9), ДНК (аденин-6)-метилтрансферазы E. coli (код доступа Swiss-Prot: P0AEE8) и ДНК (цитозин-4)-метилтрансферазы B. amyloliquefaciens (код доступа Swiss-Prot: P23941).
Чтобы скачать данные последовательности был написан конвейер:
echo 'sw:P0AED9,sw:P0AEE8,sw:P23941' | tr ',' '\n' | seqret list::/dev/stdin -outseq query_MTases.fasta -osformat2 fasta
Поиск производился программой blastp по созданной базе ORFs:
blastp -db ORFs -query query_MTases.fasta -outfmt 7 -out blastp_table.txtСсылка на файл с выдачей (таблица)
Лучшей оказалась находка с весом 50.8, её идентификатор: NZ_CP010407.1_34241Б, её координаты: [7537440 - 7538360]. Данной находке соответствовала последовательность ДНК (цитозин-4)-метилтрансферазы B. amyloliquefaciens.
По указанным выше координатам проводился поиск CDS в таблице локальных особенностей (genomic.gff), а точнее по её строкам, соответствующим CDS находки, и столбцам 4 и 5 (координаты), 7 (информация о том, является ли цепь ДНК прямой или обратной), 9 (дополнительная информация), записанным в отдельный файл командой
grep "NZ_CP010407.1" genomic.gff | grep "CDS" | cut -f4,5,7,9 > CDS.tsv
Для поиска использовала
echo -e "7537440\t7538360\t-\tMY-ORF" | cat - CDS.tsv | sort -n | grep -C 3 "MY-ORF" > neighbors.tsvСсылка на файл с выдачей
Пересечение было с CDS 7537617-7538363, в соотвествующем столбце с дополнительной информации было указано, что это ген метилтрансферазы.
Для разрешения вопроса о том, можно ли найти CDS той же находки поиском по аннотации кодирующих участков в геноме, использовала команды:
elink -target protein -db nuccore -id "NZ_CP010407.1" | efilter -query "2.1.1.113" | efetch -format fasta > search1.fasta
elink -target protein -db nuccore -id "NZ_CP010407.1" | efilter -query "2.1.1.72" | efetch -format fasta > search2.fasta
elink -target protein -db nuccore -id "NZ_CP010407.1" | efilter -query "2.1.1.37" | efetch -format fasta > search3.fasta
где 2.1.1.113, 2.1.1.72 и 2.1.1.37 - EC-коды ДНК (цитозин-4)-метилтрансферазы, ДНК (аденин-6)-метилтрансферазы, ДНК (цитозин-5)-метилтрансферазы соответственно.
С помощью третьего конвейера нашелся один белок (WP_041129408.1), как и ожидалось, ДНК (цитозин-4)-метилтрансфераза. Данный белок не соответсвует находке, полученной ранее (WP_041132502.1). При использовании первого и второго конвейера белков не нашлось.