В этом практическом задании осуществляется исследование генома и протеома бактерии Coxiella burnetii. Протеом этой бактерии имеет идентификатор UP000002671. Ранее был загружен файл с данными о белках её протеома в формате swiss с именем UP000002671.swiss.gz.
Для определения TaxID организма, применялся конвейер:
zgrep 'TaxID' UP000002671.swiss.gz | cut -d " " -f4 | sort -u | uniq -c
Файл UP000002671.swiss.gz представляет собой тот же самый загруженный файл, в котором содержится аннотация белков протеома. Все значения в этом поле были одинаковыми для каждого из белков протеома. Результат оказался следующим:
1 NCBI_TaxID=227377
1 NCBI_TaxID=227377;
Получается, что TaxID моего организма - 227377. Затем при помощи конвейера я просматривала все доступные геномные сборки для этого организма:
datasets summary genome taxon 227377 --as-json-lines | dataformat tsv genome > genome_assembly.tsv
В результате данного конвейера я получила список сборок для моего организма. Из всех я выбрала референсую версии RefSeq. AC сборки: GCF_000007765.2
Затем с использованием команды я загрузила последовательность генома и таблицу локальных характеристик.:
datasets download genome accession GCF_000007765.2 --include genome,gff3 --filename genome_data.zip
Скачанный архив ncbi_datasets.zip был разархивирован следующей командой:
unzip ncbi_datasets.zip
Чтобы выяснить, какой вариант генетического кода применяется исследуемым организмом, я загрузила таблицу с данными о таксоне из базы NCBI Taxonomy с помощью следующей команды:
efetch -db 'taxonomy' -id '227377' -format 'xml' > tax.xml
В результате выполнения программы было установлено, что организм использует таблицу №11. После этого поиск открытых рамок считывания был проведен с помощью следующей команды:
getorf -sequence GCF_000007765.2_ASM776v2_genomic.fna -outseq orfs.fasta -minsize 150 -table 11
Для того чтобы убедиться, что среди трансляций нет таких, которые короче 50 аминокислот, был применен конвейер:
infoseq -sequence orfs.fasta -only -length | sort -u -n | less
Длина всех транслированных белков составила не менее 50 аминокислот. Поэтому следующим шагом я создала белковую базу данных, основываясь на полученных открытых рамках считывания:
makeblastdb -in orfs.fasta -dbtype prot -out ORFs
Далее необходимо было получить последовательности гомологичных метилтрансфераз ( Dcm, m5C-МТаза, E.coli (код доступа в Swiss-Prot P0AED9), Dam, m6A-МТаза, E.coli (код доступа в Swiss-Prot P0AEE8), m4C-МТаза, Bacillus amyloliquefaciens (код доступа в Swiss-Prot P23941) ) без создания промежуточного файла. Для этого я воспользовалась следующим конвейером:
echo "sw:P0AED9 sw:P0AEE8 sw:P23941" | tr ' ' '\n' | seqret @/dev/stdin -outseq query_MTases.fasta -osformat fasta
Далее я с помощью последовательностей белков в качестве запроса и ORFs в качестве базы данных выполнила команд, которая выполняла поиск по сходству последовательностей:
blastp -query query_MTases.fasta -db ORFs -outfmt 7 -out blast_results.txt
Результаты:
Полная выдача blast.В результате было найдено 20 находок, среди которых лучший я посчитала по критериям минимального e-value и максимального bit-score:
subject acc.ver | s. start | s. end |
E-value | bit score | гомолог |
---|---|---|---|---|---|
NC_002971.4_9890 | 581692 | 581537 | 0.10 | 28.5 | m4C |
Далее, используя приведенный ниже конвейер, я получила файл CDS.tsv, содержащий строки и столбцы, необходимые для дальнейшей работы:
grep NC_002971.4 genomic.gff | grep CDS | cut -f4,5,7,9 > CDS.tsv
Затем были выявлены участки, которые по своим координатам наиболее близки с моей находкой:
echo -e '581537\t581692\t-\tMY-ORF' | cat - CDS.tsv | sort -n | grep -C 3 'MY-ORF' > neighbors.tsv
Были получены координаты соседних CDS:
Ссылка на файл с CDS.Получилось, что найденная мной открытая рамка пересекается с одним CDS, но координаты не совпадают
Команда для поиска CDS в аннотации генома Coxiella burnetii (2.1.1.113 – EC-код ДНК-метилтрансферазы m4C):
elink -target protein -db nuccore -id NC_002971.4 | efilter -query '2.1.1.113' | efetch -format 'fasta' > cds.fasta
Я перепробовала поискать по всем EC-кодам, но ни один из команд не выдал мне результатов. Следовательно, количество найденных белков равняется нулю.